# Need help with a research project

Page 1 of 1

## 0 Replies - 860 Views - Last Post: 10 December 2016 - 09:03 PM

### #1 hiephoang

Reputation: -1
• Posts: 23
• Joined: 18-December 15

# Need help with a research project

Posted 10 December 2016 - 09:03 PM

I'm doing research concerning Heat and Mass Transfer.

The program created in order to solve a system of 3 Ordinary Differential Equations using Runge Kutta and shooting method:

S'''+3*S*S''-2*(S'*S')+P-0.5*F=0
P''+3*10*S*P'+0.5*F'*P'+0.5(P'*P')=0
F''+3*10*10*S*F'+P''=0

Initial and boundary conditions:
x=0: S=0, S'=0, P=1, F=1
x= infinity: S'=0, P=0, F=0
Is there any one who can help me with this program? The results produced contain some errors as following:

```number of interval in the coordinate n= 300000
step in the coordinate h= 3.33333e-005
bb= 2.84057
cc= 0.777512
ee= -3.45846e-323
iter = 1 dd= 3.61808
bb= -1.#IND
cc= -1.#IND
ee= -1.#IND
iter = 2 dd= 1.#QNAN
x[i] = 0then u[i]= 0
x[i] = 0then y[i]= 0
x[i] = 0then z[i]= -1.#IND
x[i] = 0then w[i]= 1
x[i] = 0then p[i]= -1.#IND
x[i] = 0then m[i]= 1
x[i] = 0then v[i]= -1.#IND
x[i] = 3.33333e-005then u[i]= -1.#IND
x[i] = 3.33333e-005then y[i]= -1.#IND
x[i] = 3.33333e-005then z[i]= -1.#IND
x[i] = 3.33333e-005then w[i]= -1.#IND
x[i] = 3.33333e-005then p[i]= -1.#IND
x[i] = 3.33333e-005then m[i]= -1.#IND
x[i] = 3.33333e-005then v[i]= -1.#IND
x[i] = 6.66667e-005then u[i]= -1.#IND
x[i] = 6.66667e-005then y[i]= -1.#IND
x[i] = 6.66667e-005then z[i]= -1.#IND
x[i] = 6.66667e-005then w[i]= -1.#IND
x[i] = 6.66667e-005then p[i]= -1.#IND
x[i] = 6.66667e-005then m[i]= -1.#IND
x[i] = 6.66667e-005then v[i]= -1.#IND
x[i] = 0.0001then u[i]= -1.#IND
x[i] = 0.0001then y[i]= -1.#IND
x[i] = 0.0001then z[i]= -1.#IND
x[i] = 0.0001then w[i]= -1.#IND
x[i] = 0.0001then p[i]= -1.#IND
x[i] = 0.0001then m[i]= -1.#IND
x[i] = 0.0001then v[i]= -1.#IND
x[i] = 0.000133333then u[i]= -1.#IND
x[i] = 0.000133333then y[i]= -1.#IND
x[i] = 0.000133333then z[i]= -1.#IND
x[i] = 0.000133333then w[i]= -1.#IND
x[i] = 0.000133333then p[i]= -1.#IND
x[i] = 0.000133333then m[i]= -1.#IND
x[i] = 0.000133333then v[i]= -1.#IND
```

The code is following:
```#include<math.h>
double f1(double e1, double e2, double e3, double e4, double e5, double e6, double e7, double e8)
{
return e3;
}
//---------------------
double f2(double e1, double e2, double e3, double e4, double e5, double e6, double e7, double e8)
{
return e4;
}
//---------------------
double f3(double e1, double e2, double e3, double e4, double e5, double e6, double e7, double e8)
{
return  -3.0*e2*e4 + 2.0*e3*e3 - e5 + 0.5*e7;
}
//---------------------
double f4(double e1, double e2, double e3, double e4, double e5, double e6, double e7, double e8)
{
return e6;
}
double f5(double e1, double e2, double e3, double e4, double e5, double e6, double e7, double e8)
{
return -30.0*e2*e6-0.5*e8*e6-0.5*e6*e6;
}
double f6(double e1, double e2, double e3, double e4, double e5, double e6, double e7, double e8)
{
return e8;
}
double f7(double e1, double e2, double e3, double e4, double e5, double e6, double e7, double e8)
{
return -3.0*10*10*e2*e8+3.0*10*e2*e6+0.5*e8*e6+0.5*e6*e6;
}
#include<stdio.h>
#include<conio.h>
#include<iostream>
#include<fstream>
using namespace std;

//------------------------
int main(void)
{
ofstream scientificwork2("scientificwork2.txt");
const double EPS =1e-10;
int i, n = 300000;
double h=10.0/n;
scientificwork2 << "number of interval in the coordinate n= " << n << endl;
//
//
//
double *x;
x = new double[n + 1];
x[0] = 0.0;
//
for (i = 0; i<n + 1; i++)
x[i] = i*h;
//
//
scientificwork2<<"step in the coordinate h= "<<h<<endl;
//
double *u;
u = new double[n + 1];
u[0] = 0.0;
double *y;
y = new double[n + 1];
y[0] = 0.0;
double *z;
z = new double[n + 1];
z[0] = 0.7;
double *w;
w = new double[n + 1];
w[0] = 1.0;
double *p;
p = new double[n + 1];
p[0] = -0.3;
double *m;
m = new double[n + 1];
m[0] = 1.0;
double *v;
v = new double[n + 1];
v[0] = -1.25;
//
int iter;
//
double eta = 0.7, ep = -0.3, ha = -1.25;
double hz, hp, hv, dd, dz, dp, dv;
double bb, cc, ee, bbz, ccz, eez, bbp, ccp, eep, bbv, ccv, eev ;
double a1, a2, a3, a4;
double b1, b2, b3, b4;
double c1, c2, c3, c4;
double d1, d2, d3, d4;
double g1, g2, g3, g4;
double r1, r2, r3, r4;
double t1, t2, t3, t4;
iter = 0;
do
{
iter++;
u[0] = 0.0;
y[0] = 0.0;
z[0] = eta;
w[0] = 1.0;
p[0] = ep;
m[0] = 1.0;
v[0] = ha;
//
//
for (i = 0; i<n; i++)

{
a1 = f1(x[i], u[i], y[i], z[i], w[i], p[i], m[i], v[i])*h;
b1 = f2(x[i], u[i], y[i], z[i], w[i], p[i], m[i], v[i])*h;
c1 = f3(x[i], u[i], y[i], z[i], w[i], p[i], m[i], v[i])*h;
d1 = f4(x[i], u[i], y[i], z[i], w[i], p[i], m[i], v[i])*h;
g1 = f5(x[i], u[i], y[i], z[i], w[i], p[i], m[i], v[i])*h;
r1 = f6(x[i], u[i], y[i], z[i], w[i], p[i], m[i], v[i])*h;
t1 = f7(x[i], u[i], y[i], z[i], w[i], p[i], m[i], v[i])*h;
//
a2 = f1(x[i] + 0.5*h, u[i] + 0.5*a1, y[i] + 0.5*b1, z[i] + 0.5*c1, w[i] + 0.5*d1, p[i] + 0.5*g1, m[i] + 0.5*r1, v[i] + 0.5*t1)*h;
b2 = f2(x[i] + 0.5*h, u[i] + 0.5*a1, y[i] + 0.5*b1, z[i] + 0.5*c1, w[i] + 0.5*d1, p[i] + 0.5*g1, m[i] + 0.5*r1, v[i] + 0.5*t1)*h;
c2 = f3(x[i] + 0.5*h, u[i] + 0.5*a1, y[i] + 0.5*b1, z[i] + 0.5*c1, w[i] + 0.5*d1, p[i] + 0.5*g1, m[i] + 0.5*r1, v[i] + 0.5*t1)*h;
d2 = f4(x[i] + 0.5*h, u[i] + 0.5*a1, y[i] + 0.5*b1, z[i] + 0.5*c1, w[i] + 0.5*d1, p[i] + 0.5*g1, m[i] + 0.5*r1, v[i] + 0.5*t1)*h;
g2 = f5(x[i] + 0.5*h, u[i] + 0.5*a1, y[i] + 0.5*b1, z[i] + 0.5*c1, w[i] + 0.5*d1, p[i] + 0.5*g1, m[i] + 0.5*r1, v[i] + 0.5*t1)*h;
r2 = f6(x[i] + 0.5*h, u[i] + 0.5*a1, y[i] + 0.5*b1, z[i] + 0.5*c1, w[i] + 0.5*d1, p[i] + 0.5*g1, m[i] + 0.5*r1, v[i] + 0.5*t1)*h;
t2 = f7(x[i] + 0.5*h, u[i] + 0.5*a1, y[i] + 0.5*b1, z[i] + 0.5*c1, w[i] + 0.5*d1, p[i] + 0.5*g1, m[i] + 0.5*r1, v[i] + 0.5*t1)*h;
//
//‡		f1 returned	0.25000000000000000	double
//
a3 = f1(x[i] + 0.5*h, u[i] + 0.5*a2, y[i] + 0.5*b2, z[i] + 0.5*c2, w[i] + 0.5*d2, p[i] + 0.5*g2, m[i] + 0.5*r2, v[i] + 0.5*t2)*h;
b3 = f2(x[i] + 0.5*h, u[i] + 0.5*a2, y[i] + 0.5*b2, z[i] + 0.5*c2, w[i] + 0.5*d2, p[i] + 0.5*g2, m[i] + 0.5*r2, v[i] + 0.5*t2)*h;
c3 = f3(x[i] + 0.5*h, u[i] + 0.5*a2, y[i] + 0.5*b2, z[i] + 0.5*c2, w[i] + 0.5*d2, p[i] + 0.5*g2, m[i] + 0.5*r2, v[i] + 0.5*t2)*h;
d3 = f4(x[i] + 0.5*h, u[i] + 0.5*a2, y[i] + 0.5*b2, z[i] + 0.5*c2, w[i] + 0.5*d2, p[i] + 0.5*g2, m[i] + 0.5*r2, v[i] + 0.5*t2)*h;
g3 = f5(x[i] + 0.5*h, u[i] + 0.5*a2, y[i] + 0.5*b2, z[i] + 0.5*c2, w[i] + 0.5*d2, p[i] + 0.5*g2, m[i] + 0.5*r2, v[i] + 0.5*t2)*h;
r3 = f6(x[i] + 0.5*h, u[i] + 0.5*a2, y[i] + 0.5*b2, z[i] + 0.5*c2, w[i] + 0.5*d2, p[i] + 0.5*g2, m[i] + 0.5*r2, v[i] + 0.5*t2)*h;
t3 = f7(x[i] + 0.5*h, u[i] + 0.5*a2, y[i] + 0.5*b2, z[i] + 0.5*c2, w[i] + 0.5*d2, p[i] + 0.5*g2, m[i] + 0.5*r2, v[i] + 0.5*t2)*h;
//
//
a4 = f1(x[i] + h, u[i] + a3, y[i] + b3, z[i] + c3, w[i] + d3, p[i] + g3, m[i] + r3, v[i] + t3)*h;
b4 = f2(x[i] + h, u[i] + a3, y[i] + b3, z[i] + c3, w[i] + d3, p[i] + g3, m[i] + r3, v[i] + t3)*h;
c4 = f3(x[i] + h, u[i] + a3, y[i] + b3, z[i] + c3, w[i] + d3, p[i] + g3, m[i] + r3, v[i] + t3)*h;
d4 = f4(x[i] + h, u[i] + a3, y[i] + b3, z[i] + c3, w[i] + d3, p[i] + g3, m[i] + r3, v[i] + t3)*h;
g4 = f5(x[i] + h, u[i] + a3, y[i] + b3, z[i] + c3, w[i] + d3, p[i] + g3, m[i] + r3, v[i] + t3)*h;
r4 = f6(x[i] + h, u[i] + a3, y[i] + b3, z[i] + c3, w[i] + d3, p[i] + g3, m[i] + r3, v[i] + t3)*h;
t4 = f7(x[i] + h, u[i] + a3, y[i] + b3, z[i] + c3, w[i] + d3, p[i] + g3, m[i] + r3, v[i] + t3)*h;

//
u[i + 1] = u[i] + (a1 + 2.0*a2 + 2.0*a3 + a4) / 6.0;
y[i + 1] = y[i] + (b1 + 2.0*b2 + 2.0*b3 + b4) / 6.0;
z[i + 1] = z[i] + (c1 + 2.0*c2 + 2.0*c3 + c4) / 6.0;
w[i + 1] = w[i] + (d1 + 2.0*d2 + 2.0*d3 + d4) / 6.0;
p[i + 1] = p[i] + (g1 + 2.0*g2 + 2.0*g3 + g4) / 6.0;
m[i + 1] = m[i] + (r1 + 2.0*r2 + 2.0*r3 + r4) / 6.0;
v[i + 1] = v[i] + (t1 + 2.0*t2 + 2.0*t3 + t4) / 6.0;
}
//
bb = y[n] - 0.0;
cc = w[n] - 0.0;
ee = v[n] - 0.0;
scientificwork2 << "bb= " << bb << endl;
scientificwork2 << "cc= " << cc << endl;
scientificwork2 << "ee= " << ee << endl;
//
//
u[0] = 0.0;
y[0] = 0.0;
if (fabs(eta) <= 1.0) hz = sqrt(EPS);
if (fabs(eta)>1.0) hz = eta*sqrt(EPS);
z[0] = eta + hz;
w[0] = 1.0;
p[0] = ep;
m[0] = 1.0;
v[0] = ha;
//
//
//
for (i = 0; i<n; i++)

{
a1 = f1(x[i], u[i], y[i], z[i], w[i], p[i], m[i], v[i])*h;
b1 = f2(x[i], u[i], y[i], z[i], w[i], p[i], m[i], v[i])*h;
c1 = f3(x[i], u[i], y[i], z[i], w[i], p[i], m[i], v[i])*h;
d1 = f4(x[i], u[i], y[i], z[i], w[i], p[i], m[i], v[i])*h;
g1 = f5(x[i], u[i], y[i], z[i], w[i], p[i], m[i], v[i])*h;
r1 = f6(x[i], u[i], y[i], z[i], w[i], p[i], m[i], v[i])*h;
t1 = f7(x[i], u[i], y[i], z[i], w[i], p[i], m[i], v[i])*h;
//
a2 = f1(x[i] + 0.5*h, u[i] + 0.5*a1, y[i] + 0.5*b1, z[i] + 0.5*c1, w[i] + 0.5*d1, p[i] + 0.5*g1, m[i] + 0.5*r1, v[i] + 0.5*t1)*h;
b2 = f2(x[i] + 0.5*h, u[i] + 0.5*a1, y[i] + 0.5*b1, z[i] + 0.5*c1, w[i] + 0.5*d1, p[i] + 0.5*g1, m[i] + 0.5*r1, v[i] + 0.5*t1)*h;
c2 = f3(x[i] + 0.5*h, u[i] + 0.5*a1, y[i] + 0.5*b1, z[i] + 0.5*c1, w[i] + 0.5*d1, p[i] + 0.5*g1, m[i] + 0.5*r1, v[i] + 0.5*t1)*h;
d2 = f4(x[i] + 0.5*h, u[i] + 0.5*a1, y[i] + 0.5*b1, z[i] + 0.5*c1, w[i] + 0.5*d1, p[i] + 0.5*g1, m[i] + 0.5*r1, v[i] + 0.5*t1)*h;
g2 = f5(x[i] + 0.5*h, u[i] + 0.5*a1, y[i] + 0.5*b1, z[i] + 0.5*c1, w[i] + 0.5*d1, p[i] + 0.5*g1, m[i] + 0.5*r1, v[i] + 0.5*t1)*h;
r2 = f6(x[i] + 0.5*h, u[i] + 0.5*a1, y[i] + 0.5*b1, z[i] + 0.5*c1, w[i] + 0.5*d1, p[i] + 0.5*g1, m[i] + 0.5*r1, v[i] + 0.5*t1)*h;
t2 = f7(x[i] + 0.5*h, u[i] + 0.5*a1, y[i] + 0.5*b1, z[i] + 0.5*c1, w[i] + 0.5*d1, p[i] + 0.5*g1, m[i] + 0.5*r1, v[i] + 0.5*t1)*h;
//
//‡		f1 returned	0.25000000000000000	double
//
a3 = f1(x[i] + 0.5*h, u[i] + 0.5*a2, y[i] + 0.5*b2, z[i] + 0.5*c2, w[i] + 0.5*d2, p[i] + 0.5*g2, m[i] + 0.5*r2, v[i] + 0.5*t2)*h;
b3 = f2(x[i] + 0.5*h, u[i] + 0.5*a2, y[i] + 0.5*b2, z[i] + 0.5*c2, w[i] + 0.5*d2, p[i] + 0.5*g2, m[i] + 0.5*r2, v[i] + 0.5*t2)*h;
c3 = f3(x[i] + 0.5*h, u[i] + 0.5*a2, y[i] + 0.5*b2, z[i] + 0.5*c2, w[i] + 0.5*d2, p[i] + 0.5*g2, m[i] + 0.5*r2, v[i] + 0.5*t2)*h;
d3 = f4(x[i] + 0.5*h, u[i] + 0.5*a2, y[i] + 0.5*b2, z[i] + 0.5*c2, w[i] + 0.5*d2, p[i] + 0.5*g2, m[i] + 0.5*r2, v[i] + 0.5*t2)*h;
g3 = f5(x[i] + 0.5*h, u[i] + 0.5*a2, y[i] + 0.5*b2, z[i] + 0.5*c2, w[i] + 0.5*d2, p[i] + 0.5*g2, m[i] + 0.5*r2, v[i] + 0.5*t2)*h;
r3 = f6(x[i] + 0.5*h, u[i] + 0.5*a2, y[i] + 0.5*b2, z[i] + 0.5*c2, w[i] + 0.5*d2, p[i] + 0.5*g2, m[i] + 0.5*r2, v[i] + 0.5*t2)*h;
t3 = f7(x[i] + 0.5*h, u[i] + 0.5*a2, y[i] + 0.5*b2, z[i] + 0.5*c2, w[i] + 0.5*d2, p[i] + 0.5*g2, m[i] + 0.5*r2, v[i] + 0.5*t2)*h;
//
//
a4 = f1(x[i] + h, u[i] + a3, y[i] + b3, z[i] + c3, w[i] + d3, p[i] + g3, m[i] + r3, v[i] + t3)*h;
b4 = f2(x[i] + h, u[i] + a3, y[i] + b3, z[i] + c3, w[i] + d3, p[i] + g3, m[i] + r3, v[i] + t3)*h;
c4 = f3(x[i] + h, u[i] + a3, y[i] + b3, z[i] + c3, w[i] + d3, p[i] + g3, m[i] + r3, v[i] + t3)*h;
d4 = f4(x[i] + h, u[i] + a3, y[i] + b3, z[i] + c3, w[i] + d3, p[i] + g3, m[i] + r3, v[i] + t3)*h;
g4 = f5(x[i] + h, u[i] + a3, y[i] + b3, z[i] + c3, w[i] + d3, p[i] + g3, m[i] + r3, v[i] + t3)*h;
r4 = f6(x[i] + h, u[i] + a3, y[i] + b3, z[i] + c3, w[i] + d3, p[i] + g3, m[i] + r3, v[i] + t3)*h;
t4 = f7(x[i] + h, u[i] + a3, y[i] + b3, z[i] + c3, w[i] + d3, p[i] + g3, m[i] + r3, v[i] + t3)*h;
//
u[i + 1] = u[i] + (a1 + 2.0*a2 + 2.0*a3 + a4) / 6.0;
y[i + 1] = y[i] + (b1 + 2.0*b2 + 2.0*b3 + b4) / 6.0;
z[i + 1] = z[i] + (c1 + 2.0*c2 + 2.0*c3 + c4) / 6.0;
w[i + 1] = w[i] + (d1 + 2.0*d2 + 2.0*d3 + d4) / 6.0;
p[i + 1] = p[i] + (g1 + 2.0*g2 + 2.0*g3 + g4) / 6.0;
m[i + 1] = m[i] + (r1 + 2.0*r2 + 2.0*r3 + r4) / 6.0;
v[i + 1] = v[i] + (t1 + 2.0*t2 + 2.0*t3 + t4) / 6.0;
}
bbz = y[n] - 0.0;
ccz = w[n] - 0.0;
eez = v[n] - 0.0;
//
//
u[0] = 0.0;
y[0] = 0.0;
z[0] = eta;
w[0] = 1.0;
if (fabs(ep) <= 1.0) hp = sqrt(EPS);
if (fabs(ep)>1.0) hp = ep*sqrt(EPS);
p[0] = ep + hp;
m[0] = 1.0;
v[0] = ha;
//
//
for (i = 0; i<n; i++)

{
a1 = f1(x[i], u[i], y[i], z[i], w[i], p[i], m[i], v[i])*h;
b1 = f2(x[i], u[i], y[i], z[i], w[i], p[i], m[i], v[i])*h;
c1 = f3(x[i], u[i], y[i], z[i], w[i], p[i], m[i], v[i])*h;
d1 = f4(x[i], u[i], y[i], z[i], w[i], p[i], m[i], v[i])*h;
g1 = f5(x[i], u[i], y[i], z[i], w[i], p[i], m[i], v[i])*h;
r1 = f6(x[i], u[i], y[i], z[i], w[i], p[i], m[i], v[i])*h;
t1 = f7(x[i], u[i], y[i], z[i], w[i], p[i], m[i], v[i])*h;
//
a2 = f1(x[i] + 0.5*h, u[i] + 0.5*a1, y[i] + 0.5*b1, z[i] + 0.5*c1, w[i] + 0.5*d1, p[i] + 0.5*g1, m[i] + 0.5*r1, v[i] + 0.5*t1)*h;
b2 = f2(x[i] + 0.5*h, u[i] + 0.5*a1, y[i] + 0.5*b1, z[i] + 0.5*c1, w[i] + 0.5*d1, p[i] + 0.5*g1, m[i] + 0.5*r1, v[i] + 0.5*t1)*h;
c2 = f3(x[i] + 0.5*h, u[i] + 0.5*a1, y[i] + 0.5*b1, z[i] + 0.5*c1, w[i] + 0.5*d1, p[i] + 0.5*g1, m[i] + 0.5*r1, v[i] + 0.5*t1)*h;
d2 = f4(x[i] + 0.5*h, u[i] + 0.5*a1, y[i] + 0.5*b1, z[i] + 0.5*c1, w[i] + 0.5*d1, p[i] + 0.5*g1, m[i] + 0.5*r1, v[i] + 0.5*t1)*h;
g2 = f5(x[i] + 0.5*h, u[i] + 0.5*a1, y[i] + 0.5*b1, z[i] + 0.5*c1, w[i] + 0.5*d1, p[i] + 0.5*g1, m[i] + 0.5*r1, v[i] + 0.5*t1)*h;
r2 = f6(x[i] + 0.5*h, u[i] + 0.5*a1, y[i] + 0.5*b1, z[i] + 0.5*c1, w[i] + 0.5*d1, p[i] + 0.5*g1, m[i] + 0.5*r1, v[i] + 0.5*t1)*h;
t2 = f7(x[i] + 0.5*h, u[i] + 0.5*a1, y[i] + 0.5*b1, z[i] + 0.5*c1, w[i] + 0.5*d1, p[i] + 0.5*g1, m[i] + 0.5*r1, v[i] + 0.5*t1)*h;
//
//‡		f1 returned	0.25000000000000000	double
//
a3 = f1(x[i] + 0.5*h, u[i] + 0.5*a2, y[i] + 0.5*b2, z[i] + 0.5*c2, w[i] + 0.5*d2, p[i] + 0.5*g2, m[i] + 0.5*r2, v[i] + 0.5*t2)*h;
b3 = f2(x[i] + 0.5*h, u[i] + 0.5*a2, y[i] + 0.5*b2, z[i] + 0.5*c2, w[i] + 0.5*d2, p[i] + 0.5*g2, m[i] + 0.5*r2, v[i] + 0.5*t2)*h;
c3 = f3(x[i] + 0.5*h, u[i] + 0.5*a2, y[i] + 0.5*b2, z[i] + 0.5*c2, w[i] + 0.5*d2, p[i] + 0.5*g2, m[i] + 0.5*r2, v[i] + 0.5*t2)*h;
d3 = f4(x[i] + 0.5*h, u[i] + 0.5*a2, y[i] + 0.5*b2, z[i] + 0.5*c2, w[i] + 0.5*d2, p[i] + 0.5*g2, m[i] + 0.5*r2, v[i] + 0.5*t2)*h;
g3 = f5(x[i] + 0.5*h, u[i] + 0.5*a2, y[i] + 0.5*b2, z[i] + 0.5*c2, w[i] + 0.5*d2, p[i] + 0.5*g2, m[i] + 0.5*r2, v[i] + 0.5*t2)*h;
r3 = f6(x[i] + 0.5*h, u[i] + 0.5*a2, y[i] + 0.5*b2, z[i] + 0.5*c2, w[i] + 0.5*d2, p[i] + 0.5*g2, m[i] + 0.5*r2, v[i] + 0.5*t2)*h;
t3 = f7(x[i] + 0.5*h, u[i] + 0.5*a2, y[i] + 0.5*b2, z[i] + 0.5*c2, w[i] + 0.5*d2, p[i] + 0.5*g2, m[i] + 0.5*r2, v[i] + 0.5*t2)*h;
//
//
a4 = f1(x[i] + h, u[i] + a3, y[i] + b3, z[i] + c3, w[i] + d3, p[i] + g3, m[i] + r3, v[i] + t3)*h;
b4 = f2(x[i] + h, u[i] + a3, y[i] + b3, z[i] + c3, w[i] + d3, p[i] + g3, m[i] + r3, v[i] + t3)*h;
c4 = f3(x[i] + h, u[i] + a3, y[i] + b3, z[i] + c3, w[i] + d3, p[i] + g3, m[i] + r3, v[i] + t3)*h;
d4 = f4(x[i] + h, u[i] + a3, y[i] + b3, z[i] + c3, w[i] + d3, p[i] + g3, m[i] + r3, v[i] + t3)*h;
g4 = f5(x[i] + h, u[i] + a3, y[i] + b3, z[i] + c3, w[i] + d3, p[i] + g3, m[i] + r3, v[i] + t3)*h;
r4 = f6(x[i] + h, u[i] + a3, y[i] + b3, z[i] + c3, w[i] + d3, p[i] + g3, m[i] + r3, v[i] + t3)*h;
t4 = f7(x[i] + h, u[i] + a3, y[i] + b3, z[i] + c3, w[i] + d3, p[i] + g3, m[i] + r3, v[i] + t3)*h;
//
u[i + 1] = u[i] + (a1 + 2.0*a2 + 2.0*a3 + a4) / 6.0;
y[i + 1] = y[i] + (b1 + 2.0*b2 + 2.0*b3 + b4) / 6.0;
z[i + 1] = z[i] + (c1 + 2.0*c2 + 2.0*c3 + c4) / 6.0;
w[i + 1] = w[i] + (d1 + 2.0*d2 + 2.0*d3 + d4) / 6.0;
p[i + 1] = p[i] + (g1 + 2.0*g2 + 2.0*g3 + g4) / 6.0;
m[i + 1] = m[i] + (r1 + 2.0*r2 + 2.0*r3 + r4) / 6.0;
v[i + 1] = v[i] + (t1 + 2.0*t2 + 2.0*t3 + t4) / 6.0;
}
//
//
bbp = y[n] - 0.0;
ccp = w[n] - 0.0;
eep = v[n] - 0.0;
//
//
u[0] = 0.0;
y[0] = 0.0;
z[0] = eta;
w[0] = 1.0;
p[0] = ep;
m[0] = 1.0;
if (fabs(ha) <= 1.0) hv = sqrt(EPS);
if (fabs(ha) > 1.0) hv = ha*sqrt(EPS);
v[0] = ha + hv;
//
//
for (i = 0; i<n; i++)

{
a1 = f1(x[i], u[i], y[i], z[i], w[i], p[i], m[i], v[i])*h;
b1 = f2(x[i], u[i], y[i], z[i], w[i], p[i], m[i], v[i])*h;
c1 = f3(x[i], u[i], y[i], z[i], w[i], p[i], m[i], v[i])*h;
d1 = f4(x[i], u[i], y[i], z[i], w[i], p[i], m[i], v[i])*h;
g1 = f5(x[i], u[i], y[i], z[i], w[i], p[i], m[i], v[i])*h;
r1 = f6(x[i], u[i], y[i], z[i], w[i], p[i], m[i], v[i])*h;
t1 = f7(x[i], u[i], y[i], z[i], w[i], p[i], m[i], v[i])*h;
//
a2 = f1(x[i] + 0.5*h, u[i] + 0.5*a1, y[i] + 0.5*b1, z[i] + 0.5*c1, w[i] + 0.5*d1, p[i] + 0.5*g1, m[i] + 0.5*r1, v[i] + 0.5*t1)*h;
b2 = f2(x[i] + 0.5*h, u[i] + 0.5*a1, y[i] + 0.5*b1, z[i] + 0.5*c1, w[i] + 0.5*d1, p[i] + 0.5*g1, m[i] + 0.5*r1, v[i] + 0.5*t1)*h;
c2 = f3(x[i] + 0.5*h, u[i] + 0.5*a1, y[i] + 0.5*b1, z[i] + 0.5*c1, w[i] + 0.5*d1, p[i] + 0.5*g1, m[i] + 0.5*r1, v[i] + 0.5*t1)*h;
d2 = f4(x[i] + 0.5*h, u[i] + 0.5*a1, y[i] + 0.5*b1, z[i] + 0.5*c1, w[i] + 0.5*d1, p[i] + 0.5*g1, m[i] + 0.5*r1, v[i] + 0.5*t1)*h;
g2 = f5(x[i] + 0.5*h, u[i] + 0.5*a1, y[i] + 0.5*b1, z[i] + 0.5*c1, w[i] + 0.5*d1, p[i] + 0.5*g1, m[i] + 0.5*r1, v[i] + 0.5*t1)*h;
r2 = f6(x[i] + 0.5*h, u[i] + 0.5*a1, y[i] + 0.5*b1, z[i] + 0.5*c1, w[i] + 0.5*d1, p[i] + 0.5*g1, m[i] + 0.5*r1, v[i] + 0.5*t1)*h;
t2 = f7(x[i] + 0.5*h, u[i] + 0.5*a1, y[i] + 0.5*b1, z[i] + 0.5*c1, w[i] + 0.5*d1, p[i] + 0.5*g1, m[i] + 0.5*r1, v[i] + 0.5*t1)*h;
//
//‡		f1 returned	0.25000000000000000	double
//
a3 = f1(x[i] + 0.5*h, u[i] + 0.5*a2, y[i] + 0.5*b2, z[i] + 0.5*c2, w[i] + 0.5*d2, p[i] + 0.5*g2, m[i] + 0.5*r2, v[i] + 0.5*t2)*h;
b3 = f2(x[i] + 0.5*h, u[i] + 0.5*a2, y[i] + 0.5*b2, z[i] + 0.5*c2, w[i] + 0.5*d2, p[i] + 0.5*g2, m[i] + 0.5*r2, v[i] + 0.5*t2)*h;
c3 = f3(x[i] + 0.5*h, u[i] + 0.5*a2, y[i] + 0.5*b2, z[i] + 0.5*c2, w[i] + 0.5*d2, p[i] + 0.5*g2, m[i] + 0.5*r2, v[i] + 0.5*t2)*h;
d3 = f4(x[i] + 0.5*h, u[i] + 0.5*a2, y[i] + 0.5*b2, z[i] + 0.5*c2, w[i] + 0.5*d2, p[i] + 0.5*g2, m[i] + 0.5*r2, v[i] + 0.5*t2)*h;
g3 = f5(x[i] + 0.5*h, u[i] + 0.5*a2, y[i] + 0.5*b2, z[i] + 0.5*c2, w[i] + 0.5*d2, p[i] + 0.5*g2, m[i] + 0.5*r2, v[i] + 0.5*t2)*h;
r3 = f6(x[i] + 0.5*h, u[i] + 0.5*a2, y[i] + 0.5*b2, z[i] + 0.5*c2, w[i] + 0.5*d2, p[i] + 0.5*g2, m[i] + 0.5*r2, v[i] + 0.5*t2)*h;
t3 = f7(x[i] + 0.5*h, u[i] + 0.5*a2, y[i] + 0.5*b2, z[i] + 0.5*c2, w[i] + 0.5*d2, p[i] + 0.5*g2, m[i] + 0.5*r2, v[i] + 0.5*t2)*h;
//
//

a4 = f1(x[i] + h, u[i] + a3, y[i] + b3, z[i] + c3, w[i] + d3, p[i] + g3, m[i] + r3, v[i] + t3)*h;
b4 = f2(x[i] + h, u[i] + a3, y[i] + b3, z[i] + c3, w[i] + d3, p[i] + g3, m[i] + r3, v[i] + t3)*h;
c4 = f3(x[i] + h, u[i] + a3, y[i] + b3, z[i] + c3, w[i] + d3, p[i] + g3, m[i] + r3, v[i] + t3)*h;
d4 = f4(x[i] + h, u[i] + a3, y[i] + b3, z[i] + c3, w[i] + d3, p[i] + g3, m[i] + r3, v[i] + t3)*h;
g4 = f5(x[i] + h, u[i] + a3, y[i] + b3, z[i] + c3, w[i] + d3, p[i] + g3, m[i] + r3, v[i] + t3)*h;
r4 = f6(x[i] + h, u[i] + a3, y[i] + b3, z[i] + c3, w[i] + d3, p[i] + g3, m[i] + r3, v[i] + t3)*h;
t4 = f7(x[i] + h, u[i] + a3, y[i] + b3, z[i] + c3, w[i] + d3, p[i] + g3, m[i] + r3, v[i] + t3)*h;
//
u[i + 1] = u[i] + (a1 + 2.0*a2 + 2.0*a3 + a4) / 6.0;
y[i + 1] = y[i] + (b1 + 2.0*b2 + 2.0*b3 + b4) / 6.0;
z[i + 1] = z[i] + (c1 + 2.0*c2 + 2.0*c3 + c4) / 6.0;
w[i + 1] = w[i] + (d1 + 2.0*d2 + 2.0*d3 + d4) / 6.0;
p[i + 1] = p[i] + (g1 + 2.0*g2 + 2.0*g3 + g4) / 6.0;
m[i + 1] = m[i] + (r1 + 2.0*r2 + 2.0*r3 + r4) / 6.0;
v[i + 1] = v[i] + (t1 + 2.0*t2 + 2.0*t3 + t4) / 6.0;
}
//
//
bbv = y[n] - 0.0;
ccv = w[n] - 0.0;
eev = v[n] - 0.0;
//
//
dz = (cc*(bbp - bb) / hp - bb*(ccp - cc) / hp) / ((bbz - bb) / hz*(ccp - cc) / hp - (bbp - bb) / hp*(ccz - cc) / hz);
dp = (cc*(bbz - bb) / hz - bb*(ccz - cc) / hz) / ((ccz - cc) / hz*(bbp - bb) / hp - (bbz - bb) / hz*(ccp - cc) / hp);
dv = (cc*(bbz - bb) / hz - bb*(ccz - cc) / hz) / ((ccz - cc) / hz*(bbv - bb) / hv - (bbz - bb) / hz*(ccv - cc) / hv);
//
//
eta = eta + dz;
ep = ep + dp;
ha = ha + dv;
//
//
dd = fabs(bb) + fabs(cc) + fabs(ee);
//
scientificwork2 << "iter = " << iter << " dd= " << dd << endl;
//
} while (dd>EPS);
//
//
//
u[0] = 0.0;
y[0] = 0.0;
z[0] = eta;
w[0] = 1.0;
p[0] = ep;
m[0] = 1.0;
v[0] = ha;
//
//
for (i = 0; i<n; i++)

{
a1 = f1(x[i], u[i], y[i], z[i], w[i], p[i], m[i], v[i])*h;
b1 = f2(x[i], u[i], y[i], z[i], w[i], p[i], m[i], v[i])*h;
c1 = f3(x[i], u[i], y[i], z[i], w[i], p[i], m[i], v[i])*h;
d1 = f4(x[i], u[i], y[i], z[i], w[i], p[i], m[i], v[i])*h;
g1 = f5(x[i], u[i], y[i], z[i], w[i], p[i], m[i], v[i])*h;
r1 = f6(x[i], u[i], y[i], z[i], w[i], p[i], m[i], v[i])*h;
t1 = f7(x[i], u[i], y[i], z[i], w[i], p[i], m[i], v[i])*h;
//
a2 = f1(x[i] + 0.5*h, u[i] + 0.5*a1, y[i] + 0.5*b1, z[i] + 0.5*c1, w[i] + 0.5*d1, p[i] + 0.5*g1, m[i] + 0.5*r1, v[i] + 0.5*t1)*h;
b2 = f2(x[i] + 0.5*h, u[i] + 0.5*a1, y[i] + 0.5*b1, z[i] + 0.5*c1, w[i] + 0.5*d1, p[i] + 0.5*g1, m[i] + 0.5*r1, v[i] + 0.5*t1)*h;
c2 = f3(x[i] + 0.5*h, u[i] + 0.5*a1, y[i] + 0.5*b1, z[i] + 0.5*c1, w[i] + 0.5*d1, p[i] + 0.5*g1, m[i] + 0.5*r1, v[i] + 0.5*t1)*h;
d2 = f4(x[i] + 0.5*h, u[i] + 0.5*a1, y[i] + 0.5*b1, z[i] + 0.5*c1, w[i] + 0.5*d1, p[i] + 0.5*g1, m[i] + 0.5*r1, v[i] + 0.5*t1)*h;
g2 = f5(x[i] + 0.5*h, u[i] + 0.5*a1, y[i] + 0.5*b1, z[i] + 0.5*c1, w[i] + 0.5*d1, p[i] + 0.5*g1, m[i] + 0.5*r1, v[i] + 0.5*t1)*h;
r2 = f6(x[i] + 0.5*h, u[i] + 0.5*a1, y[i] + 0.5*b1, z[i] + 0.5*c1, w[i] + 0.5*d1, p[i] + 0.5*g1, m[i] + 0.5*r1, v[i] + 0.5*t1)*h;
t2 = f7(x[i] + 0.5*h, u[i] + 0.5*a1, y[i] + 0.5*b1, z[i] + 0.5*c1, w[i] + 0.5*d1, p[i] + 0.5*g1, m[i] + 0.5*r1, v[i] + 0.5*t1)*h;
//
//‡		f1 returned	0.25000000000000000	double
//
a3 = f1(x[i] + 0.5*h, u[i] + 0.5*a2, y[i] + 0.5*b2, z[i] + 0.5*c2, w[i] + 0.5*d2, p[i] + 0.5*g2, m[i] + 0.5*r2, v[i] + 0.5*t2)*h;
b3 = f2(x[i] + 0.5*h, u[i] + 0.5*a2, y[i] + 0.5*b2, z[i] + 0.5*c2, w[i] + 0.5*d2, p[i] + 0.5*g2, m[i] + 0.5*r2, v[i] + 0.5*t2)*h;
c3 = f3(x[i] + 0.5*h, u[i] + 0.5*a2, y[i] + 0.5*b2, z[i] + 0.5*c2, w[i] + 0.5*d2, p[i] + 0.5*g2, m[i] + 0.5*r2, v[i] + 0.5*t2)*h;
d3 = f4(x[i] + 0.5*h, u[i] + 0.5*a2, y[i] + 0.5*b2, z[i] + 0.5*c2, w[i] + 0.5*d2, p[i] + 0.5*g2, m[i] + 0.5*r2, v[i] + 0.5*t2)*h;
g3 = f5(x[i] + 0.5*h, u[i] + 0.5*a2, y[i] + 0.5*b2, z[i] + 0.5*c2, w[i] + 0.5*d2, p[i] + 0.5*g2, m[i] + 0.5*r2, v[i] + 0.5*t2)*h;
r3 = f6(x[i] + 0.5*h, u[i] + 0.5*a2, y[i] + 0.5*b2, z[i] + 0.5*c2, w[i] + 0.5*d2, p[i] + 0.5*g2, m[i] + 0.5*r2, v[i] + 0.5*t2)*h;
t3 = f7(x[i] + 0.5*h, u[i] + 0.5*a2, y[i] + 0.5*b2, z[i] + 0.5*c2, w[i] + 0.5*d2, p[i] + 0.5*g2, m[i] + 0.5*r2, v[i] + 0.5*t2)*h;
//
//
a4 = f1(x[i] + h, u[i] + a3, y[i] + b3, z[i] + c3, w[i] + d3, p[i] + g3, m[i] + r3, v[i] + t3)*h;
b4 = f2(x[i] + h, u[i] + a3, y[i] + b3, z[i] + c3, w[i] + d3, p[i] + g3, m[i] + r3, v[i] + t3)*h;
c4 = f3(x[i] + h, u[i] + a3, y[i] + b3, z[i] + c3, w[i] + d3, p[i] + g3, m[i] + r3, v[i] + t3)*h;
d4 = f4(x[i] + h, u[i] + a3, y[i] + b3, z[i] + c3, w[i] + d3, p[i] + g3, m[i] + r3, v[i] + t3)*h;
g4 = f5(x[i] + h, u[i] + a3, y[i] + b3, z[i] + c3, w[i] + d3, p[i] + g3, m[i] + r3, v[i] + t3)*h;
r4 = f6(x[i] + h, u[i] + a3, y[i] + b3, z[i] + c3, w[i] + d3, p[i] + g3, m[i] + r3, v[i] + t3)*h;
t4 = f7(x[i] + h, u[i] + a3, y[i] + b3, z[i] + c3, w[i] + d3, p[i] + g3, m[i] + r3, v[i] + t3)*h;

//
u[i + 1] = u[i] + (a1 + 2.0*a2 + 2.0*a3 + a4) / 6.0;
y[i + 1] = y[i] + (b1 + 2.0*b2 + 2.0*b3 + b4) / 6.0;
z[i + 1] = z[i] + (c1 + 2.0*c2 + 2.0*c3 + c4) / 6.0;
w[i + 1] = w[i] + (d1 + 2.0*d2 + 2.0*d3 + d4) / 6.0;
p[i + 1] = p[i] + (g1 + 2.0*g2 + 2.0*g3 + g4) / 6.0;
m[i + 1] = m[i] + (r1 + 2.0*r2 + 2.0*r3 + r4) / 6.0;
v[i + 1] = v[i] + (t1 + 2.0*t2 + 2.0*t3 + t4) / 6.0;
}
//
//
for (i = 0; i<n + 1; i++)
{
scientificwork2<<"x[i] = "<<x[i]<<"then u[i]= "<<u[i]<<endl;
scientificwork2 << "x[i] = " << x[i] << "then y[i]= " << y[i] << endl;
scientificwork2 << "x[i] = " << x[i] << "then z[i]= " << z[i] << endl;
scientificwork2 << "x[i] = " << x[i] << "then w[i]= " << w[i] << endl;
scientificwork2 << "x[i] = " << x[i] << "then p[i]= " << p[i] << endl;
scientificwork2 << "x[i] = " << x[i] << "then m[i]= " << m[i] << endl;
scientificwork2 << "x[i] = " << x[i] << "then v[i]= " << v[i] << endl;

}
delete[]x;
delete[]u;
delete[]y;
delete[]z;
delete[]w;
delete[]p;
delete[]m;
delete[]v;
//
scientificwork2.close();
system("pause");
}

```

Thanks
Best

Is This A Good Question/Topic? 0

Page 1 of 1

 .related ul { list-style-type: circle; font-size: 12px; font-weight: bold; } .related li { margin-bottom: 5px; background-position: left 7px !important; margin-left: -35px; } .related h2 { font-size: 18px; font-weight: bold; } .related a { color: blue; }