home
***
CD-ROM
|
disk
|
FTP
|
other
***
search
/
Club Amiga de Montreal - CAM
/
CAM_CD_1.iso
/
files
/
509.lha
/
ATS_v2.0
/
ats2.cp
< prev
next >
Wrap
Text File
|
1991-05-06
|
5KB
|
121 lines
// ret.terms[i]=complex(0.0,0.0);
//
// C++ Routines to manipulate Taylor series for solving a differential
// equation.
//
// By Dennis J. Darland
// 4/3/91
// Public Domain
// Version 1.6
// Many enhancements expected.
// PROVIDED WITH NO WARRANTY.
//
#include "ats.h"
complex radius_cp(t_s & y,complex & order)
{
register int m;
complex rm0,rm1,rm2,rm3,rm4,rcs,c;
complex t1,t2,t3,t4;
m = y.last_term - 1;
rm0 = y.terms[m]/y.terms[m-1];
cout << "rm0 = " << rm0 << "\n";
rm1 = y.terms[m-1]/y.terms[m-2];
cout << "rm1 = " << rm1 << "\n";
rm2 = y.terms[m-2]/y.terms[m-3];
cout << "rm2 = " << rm2 << "\n";
rm3 = y.terms[m-3]/y.terms[m-4];
cout << "rm3 = " << rm3 << "\n";
rm4 = y.terms[m-4]/y.terms[m-5];
cout << "rm4 = " << rm4 << "\n";
t1 = complex(12.0,0.0)*rm2*rm1*rm1*rm4;
t1 = t1 -complex(6.0,0.0)*rm3*rm3*rm2*rm4;
t1 = t1 -complex(10.0,0.0)*rm1*rm1*rm2*rm3;
t1 = t1 +complex(8.0,0.0)*rm2*rm2*rm3*rm4;
t1 = t1 +complex(10.0,0.0)*rm2*rm2*rm1*rm3;
t1 = t1 +complex(18.0,0.0)*rm1*rm3*rm3*rm4;
t1 = t1 +complex(25.0,0.0)*complex(double(m),0.0)*rm1*rm3*rm2*rm4;
t1 = t1 -complex(5.0,0.0)*complex(double(m),0.0)*rm0*rm2*rm1*rm3;
t1 = t1 +complex(12.0,0.0)*complex(double(m),0.0)*rm1*rm1*rm2*rm3;
t1 = t1 -complex(7.0,0.0)*complex(double(m),0.0)*rm2*rm2*rm1*rm3;
t1 = t1 -complex(32.0,0.0)*rm2*rm1*rm3*rm4;
t1 = t1 -complex(8.0,0.0)*complex(double(m),0.0)*rm2*rm2*rm3*rm4;
t1 = t1 +complex(5.0,0.0)*complex(double(m),0.0)*rm3*rm3*rm2*rm4;
t1 = t1 -complex(12.0,0.0)*complex(double(m),0.0)*rm3*rm3*rm4*rm1;
t1 = t1 -complex(15.0,0.0)*complex(double(m),0.0)*rm1*rm1*rm2*rm4;
t1 = t1 +complex(8.0,0.0)*complex(double(m),0.0)*rm0*rm2*rm1*rm4;
t1 = t1 -complex(3.0,0.0)*complex(double(m),0.0)*rm0*rm1*rm3*rm4;
t1 = t1 +rm1*rm3*rm4*rm0*complex(double(m),0.0)*complex(double(m),0.0);
t1 = t1 +rm2*rm1*rm3*rm0*complex(double(m),0.0)*complex(double(m),0.0);
t1 = t1 -complex(2.0,0.0)*rm2*rm1*rm4*rm0*complex(double(m),0.0)*complex(double(m),0.0);
t1 = t1 -complex(2.0,0.0)*rm1*rm1*rm2*rm3*complex(double(m),0.0)*complex(double(m),0.0);
t1 = t1 +complex(3.0,0.0)*rm2*rm1*rm1*rm4*complex(double(m),0.0)*complex(double(m),0.0);
t1 = t1 +rm2*rm2*rm1*rm3*complex(double(m),0.0)*complex(double(m),0.0);
t1 = t1 -complex(5.0,0.0)*rm2*rm1*rm3*rm4*complex(double(m),0.0)*complex(double(m),0.0);
t1 = t1 +complex(2.0,0.0)*rm1*rm3*rm3*rm4*complex(double(m),0.0)*complex(double(m),0.0);
t1 = t1 +complex(2.0,0.0)*rm2*rm2*rm3*rm4*complex(double(m),0.0)*complex(double(m),0.0);
t1 = t1 -rm3*rm3*rm2*rm4*complex(double(m),0.0)*complex(double(m),0.0);
t2 = -complex(5.0,0.0)*complex(double(m),0.0)*rm1*rm3*rm2*rm4;
t2 = t2 +complex(3.0,0.0)*complex(double(m),0.0)*rm1*rm1*rm2*rm4;
t2 = t2 +complex(10.0,0.0)*rm2*rm1*rm3*rm4;
t2 = t2 -complex(3.0,0.0)*rm2*rm1*rm1*rm4;
t2 = t2 +complex(2.0,0.0)*complex(double(m),0.0)*rm2*rm2*rm3*rm4;
t2 = t2 -complex(4.0,0.0)*rm2*rm2*rm3*rm4;
t2 = t2 +complex(2.0,0.0)*complex(double(m),0.0)*rm3*rm3*rm4*rm1;
t2 = t2 -complex(double(m),0.0)*rm3*rm3*rm2*rm4;
t2 = t2 -complex(6.0,0.0)*rm1*rm3*rm3*rm4;
t2 = t2 +complex(3.0,0.0)*rm3*rm3*rm2*rm4;
t2 = t2 -complex(2.0,0.0)*complex(double(m),0.0)*rm0*rm2*rm1*rm4;
t2 = t2 +complex(double(m),0.0)*rm0*rm1*rm3*rm4;
t2 = t2 +complex(double(m),0.0)*rm0*rm2*rm1*rm3;
t2 = t2 -complex(2.0,0.0)*complex(double(m),0.0)*rm1*rm1*rm2*rm3;
t2 = t2 +complex(2.0,0.0)*rm1*rm1*rm2*rm3;
t2 = t2 +complex(double(m),0.0)*rm2*rm2*rm1*rm3;
t2 = t2 -complex(2.0,0.0)*rm2*rm2*rm1*rm3;
order = -t1/(complex(2.0,0.0)*t2);
t3 = rm3*rm4;
t3 = t3 -complex(4.0,0.0)*rm2*rm4;
t3 = t3 +rm2*rm1;
t3 = t3 -complex(4.0,0.0)*rm1*rm3;
t3 = t3 +complex(3.0,0.0)*rm1*rm4;
t3 = t3 +complex(3.0,0.0)*rm2*rm3;
t4 = -complex(5.0,0.0)*complex(double(m),0.0)*rm1*rm3*rm2*rm4;
t4 = t4 +complex(3.0,0.0)*complex(double(m),0.0)*rm1*rm1*rm2*rm4;
t4 = t4 +complex(10.0,0.0)*rm2*rm1*rm3*rm4;
t4 = t4 -complex(3.0,0.0)*rm2*rm1*rm1*rm4;
t4 = t4 +complex(2.0,0.0)*complex(double(m),0.0)*rm2*rm2*rm3*rm4;
t4 = t4 -complex(4.0,0.0)*rm2*rm2*rm3*rm4;
t4 = t4 +complex(2.0,0.0)*complex(double(m),0.0)*rm3*rm3*rm4*rm1;
t4 = t4 -complex(double(m),0.0)*rm3*rm3*rm2*rm4;
t4 = t4 -complex(6.0,0.0)*rm1*rm3*rm3*rm4;
t4 = t4 +complex(3.0,0.0)*rm3*rm3*rm2*rm4;
t4 = t4 -complex(2.0,0.0)*complex(double(m),0.0)*rm0*rm2*rm1*rm4;
t4 = t4 +complex(double(m),0.0)*rm0*rm1*rm3*rm4;
t4 = t4 +complex(double(m),0.0)*rm0*rm2*rm1*rm3;
t4 = t4 -complex(2.0,0.0)*complex(double(m),0.0)*rm1*rm1*rm2*rm3;
t4 = t4 +complex(2.0,0.0)*rm1*rm1*rm2*rm3;
t4 = t4 +complex(double(m),0.0)*rm2*rm2*rm1*rm3;
t4 = t4 -complex(2.0,0.0)*rm2*rm2*rm1*rm3;
rcs = -t3/t4;
c = my_sqrt(rcs)*y.h;
return(c);
}
complex radius_xx(t_s & y,complex & order)
{
register int m;
complex rm0,rm1,hdrc,rcs;
m = y.last_term - 1;
rm0 = y.terms[m]/y.terms[m-1];
cout << "rm0 = " << rm0 << "\n";
rm1 = y.terms[m-1]/y.terms[m-2];
cout << "rm1 = " << rm1 << "\n";
hdrc = complex(double(m),0.0)*rm0-complex(double(m-1),0.0)*rm1;
cout << "hdrc = " << hdrc << "\n";
rcs = y.h/hdrc;
cout << "rcs = " << rcs << "\n";
order = complex(double(m),0.0)*rm0*rcs/y.h-(complex(double(m-1),0.0));
cout << "order = " << order << "\n";
return(rcs);
}