We present a program to solve by a numerical scheme a system of first-order linear differential equations of the form
A_{2}y'' + A_{1}y' + A_{0}y = fwhere the matrices A_{0} and A_{1} and A_{2} are
| 2 -1 | | -2 0 | | 2 0 | A0 = | |, A1 = | |, A2 = | | | -1 2 | | 0 1 | | 0 1 |and
| (2-2π^{2})sin(πt) - 2πcos(πt) - cos(2πt) | f = | |. | (2-4π^{2})cos(2πt) - 2πsin(2πt) - sin(πt) |With adapted initial conditions, the solution is given by
y_{1}(t) = sin(πt), y_{2} = cos(2πt)
ODESolver ode(NEWMARK,theTimeStep,theFinalTime,2); |
DMatrix<double> A0(2,2), A1(2,2), A2(2,2); ode.setMatrices(A0,A1,A2); Vect<double> u(2), v(2); u(1) = 0; u(2) = 1; v(1) = OFELI_PI; v(2) = 0; ode.setInitial(u,v); ode.setRHS("(2-2*pi^2)*sin(pi*t)-2*pi*cos(pi*t)-cos(2*pi*t)"); ode.setRHS("(2-4*pi^2)*cos(2*pi*t)-2*pi*sin(2*pi*t)-sin(pi*t)"); |
TimeLoop { A0(1,1) = 2.; A0(1,2) = -1.; A0(2,1) = -1.; A0(2,2) = 2.; A1(1,1) = -2.; A1(1,2) = 0.; A1(2,1) = 0.; A1(2,2) = 1.; A2(1,1) = 2.; A2(1,2) = 0.; A2(2,1) = 0.; A2(2,2) = 1.; ode.runOneTimeStep(); } |
cout << ode << endl; cout << "Error: " << fabs(u(1)-sin(OFELI_PI*theFinalTime)) << ", " << fabs(u(2)-cos(2*OFELI_PI*theFinalTime)) << endl; return 0; } |