Задание №2 Решение систем обыкновенных дифференциальных уравнений
#include <stdio.h> #include <math.h> const int n = 2; // укажите количество уравнений в вашем варианте
double f1 (double y[], double time) { return y[1]; // введите первую правую часть } double f2 (double y[], double time) { return (1.0-y[0]*y[0])*y[1]-y[0]; // введите вторую правую часть } //double f3 (double x[], double time) { // return; // введите третью правую часть //}
double (* p[n]) (double [], double);
void main () { const int m = 4; FILE * fil; double y[n] = {2.0,0.0 }; // укажите начальные значения параметров p[0] = f1; p[1] = f2; //p[2] = f3; double t = 0.0, tmax = 0.1, tau = 0.01; double r[m][n], yy[n]; int i; fopen_s (&fil, "zadanie2.dat", "w");
// Метод Эйлера printf (" Metod Ejler "); fprintf (fil, " Метод Эйлера "); t = 0.0; y[0] =2.0; y[1] = 0.0; // y[2] = 0.0; // укажите начальные значения параметров printf ("time = %1.2f\n ", t); fprintf (fil, "time = %1.2f\n ", t); for (i=0; i<n; i++) { printf ("y[%d] = %2.8f ", i, y[i]); fprintf (fil, "y[%d] = %2.8f ", i, y[i]); } printf ("\n"); fprintf (fil, "\n"); do { for (i=0; i<n; i++) yy[i] = y[i] + tau*(p[i](y,t)); for (i=0; i<n; i++) y[i] = yy[i]; t += tau; printf ("time = %1.2f\n ", t); fprintf (fil, "time = %1.2f\n ", t); for (i=0; i<n; i++) { printf ("y[%d] = %2.8f ", i, y[i]); fprintf (fil, "y[%d] = %2.8f ", i, y[i]); } printf ("\n"); fprintf (fil, "\n"); } while (t <= tmax); // Метод Рунге-Кутта 2 printf ("\n Metod Runge-Kutta 2 "); fprintf (fil, "\n Метод Рунге-Кутты 2 "); t = 0.0; y[0] = 2.0; y[1] =0.0; //y[2] = 0.0; // укажите начальные значения параметров double ff[n]; printf ("time = %1.2f\n ", t); fprintf (fil, "time = %1.2f\n ", t); for (i=0; i<n; i++) { printf ("y[%d] = %2.8f ", i, y[i]); fprintf (fil, "y[%d] = %2.8f ", i, y[i]); } printf ("\n"); fprintf (fil, "\n"); do { for (i=0; i<n; i++) yy[i] = y[i]+0.5*tau*(p[i](y,t)); for (i=0; i<n; i++) ff[i] = (p[i](yy,t+tau*0.5)); for (i=0; i<n; i++) y[i] += tau*ff[i]; t += tau; printf ("time = %1.2f\n ", t); fprintf (fil, "time = %1.2f\n ", t); for (i=0; i<n; i++) { printf ("y[%d] = %2.8f ", i, y[i]); fprintf (fil, "y[%d] = %2.8f ", i, y[i]); } printf ("\n"); fprintf (fil, "\n"); } while (t <= tmax); // Метод Рунге-Кутта 4 printf ("\n Metod Runge-Kutta 4 "); fprintf (fil, "\n Метод Рунге-Кутты 4 "); t = 0.0; y[0] =2.0; y[1] =0.0; //y[2] = 0.0; // укажите начальные значения параметров printf ("time = %1.2f\n ", t); fprintf (fil, "time = %1.2f\n ", t); for (i=0; i<n; i++) { printf ("y[%d] = %2.8f ", i, y[i]); fprintf (fil, "y[%d] = %2.8f ", i, y[i]); } printf ("\n"); fprintf (fil, "\n"); do { for (i=0; i<n; i++) r[0][i] = tau*(p[i](y,t)); for (i=0; i<n; i++) yy[i] = y[i]+0.5*r[0][i]; for (i=0; i<n; i++) r[1][i] = tau*(p[i](yy,t+0.5*tau)); for (i=0; i<n; i++) yy[i] = y[i]+0.5*r[1][i]; for (i=0; i<n; i++) r[2][i] = tau*(p[i](yy,t+0.5*tau)); for (i=0; i<n; i++) yy[i] = y[i]+r[2][i]; for (i=0; i<n; i++) r[3][i] = tau*(p[i](yy,t+tau)); for (i=0; i<n; i++) y[i] += (r[0][i]+2.0*r[1][i]+2.0*r[2][i]+r[3][i])/6.0; t += tau; printf ("time = %1.2f\n", t); fprintf (fil, "time = %1.2f\n ", t); for (i=0; i<n; i++) { printf ("y[%d] = %2.10f ", i, y[i]); fprintf (fil, "y[%d] = %2.8f ", i, y[i]); } printf ("\n"); fprintf (fil, "\n"); } while (t <= tmax); fclose (fil); }
решение
Метод Эйлера time = 0.00 y[0] = 2.00000000 y[1] = 0.00000000 time = 0.01 y[0] = 2.00000000 y[1] = -0.02000000 time = 0.02 y[0] = 1.99980000 y[1] = -0.03940000 time = 0.03 y[0] = 1.99940600 y[1] = -0.05821632 time = 0.04 y[0] = 1.99882384 y[1] = -0.07646527 time = 0.05 y[0] = 1.99805918 y[1] = -0.09416315 time = 0.06 y[0] = 1.99711755 y[1] = -0.11132615 time = 0.07 y[0] = 1.99600429 y[1] = -0.12797037 time = 0.08 y[0] = 1.99472459 y[1] = -0.14411173 time = 0.09 y[0] = 1.99328347 y[1] = -0.15976600 time = 0.10 y[0] = 1.99168581 y[1] = -0.17494870 time = 0.11 y[0] = 1.98993632 y[1] = -0.18967516
Метод Рунге-Кутты 2 time = 0.00 y[0] = 2.00000000 y[1] = 0.00000000 time = 0.01 y[0] = 1.99990000 y[1] = -0.01970000 time = 0.02 y[0] = 1.99960596 y[1] = -0.03881613 time = 0.03 y[0] = 1.99912364 y[1] = -0.05736441 time = 0.04 y[0] = 1.99845863 y[1] = -0.07536086 time = 0.05 y[0] = 1.99761638 y[1] = -0.09282138 time = 0.06 y[0] = 1.99660217 y[1] = -0.10976180 time = 0.07 y[0] = 1.99542111 y[1] = -0.12619779 time = 0.08 y[0] = 1.99407817 y[1] = -0.14214484 time = 0.09 y[0] = 1.99257817 y[1] = -0.15761824 time = 0.10 y[0] = 1.99092577 y[1] = -0.17263309 time = 0.11 y[0] = 1.98912548 y[1] = -0.18720422
Метод Рунге-Кутты 4 time = 0.00 y[0] = 2.00000000 y[1] = 0.00000000 time = 0.01 y[0] = 1.99990099 y[1] = -0.01970267 time = 0.02 y[0] = 1.99960789 y[1] = -0.03882136 time = 0.03 y[0] = 1.99912646 y[1] = -0.05737210 time = 0.04 y[0] = 1.99846229 y[1] = -0.07537089 time = 0.05 y[0] = 1.99762083 y[1] = -0.09283365 time = 0.06 y[0] = 1.99660735 y[1] = -0.10977620 time = 0.07 y[0] = 1.99542699 y[1] = -0.12621419 time = 0.08 y[0] = 1.99408470 y[1] = -0.14216314 time = 0.09 y[0] = 1.99258530 y[1] = -0.15763834 time = 0.10 y[0] = 1.99093346 y[1] = -0.17265487 time = 0.11 y[0] = 1.98913368 y[1] = -0.18722757
Задпние №3 решение уравнений в частных производных #include <stdio.h> #include <math.h>
const double L = 200.0, h = 1.0, tau = 0.1, u0 = 22.0; const double tn = 0.0, tmax = 10000.0, a = 1.0;
double g0 (double x) { return sin(x+0.480); // Введите функцию g0(x) } double f0 (double time) { return 0.46180;// Введите функцию f0(t) } double f1 (double time) { return 3.0*time+0.8820; // Введите функцию f1(t) }
int main () { FILE * fil; // блок подготовки int n = (int) (L / h) + 1; double x = 0.0; double * u = new double [n]; double time = tn; fopen_s (&fil, "teplo.dat", "w"); // начальные значения температур u[0] = f0 (time);; //fi1(time); for (int i=1; i<n-1; i++) { u[i] = g0 (x); x += h; } u[n-1] = f1 (time); //fi2(time);
// основной расчет double upr, uu, R = a*a*tau/(h*h); do { upr = u[0]; for (int i=1; i<n-1; i++) { uu = u[i]; u[i] += R*(upr-2.0*u[i]+u[i+1]); upr = uu; } u[0] = f0(time+tau); u[n-1] = f1(time+tau); time += tau; } while (time <= tmax); // вывод вычесленных значений температур fprintf (fil, "\n вывод вычесленных значений температур\n"); for (int i=0; i<n; i++) { fprintf (fil, "u[%d]=%2.8f\t", i, u[i]); if (!((i+1)%3)) fprintf (fil, "\n"); } fprintf (fil, "\n"); printf ("End calculus\n"); delete [] u; fclose (fil); return 0; }
Результаты вычислений вывод вычесленных значений температур u[0]=0.46180000 u[1]=60.77452300 u[2]=121.09970379 u[3]=181.44980078 u[4]=241.83727365 u[5]=302.27458393 u[6]=362.77419567 u[7]=423.34857605 u[8]=484.01019599 u[9]=544.77153080 u[10]=605.64506077 u[11]=666.64327182 u[12]=727.77865608 u[13]=789.06371257 u[14]=850.51094776 u[15]=912.13287620 u[16]=973.94202114 u[17]=1035.95091514 u[18]=1098.17210068 u[19]=1160.61813073 u[20]=1223.30156941 u[21]=1286.23499256 u[22]=1349.43098833 u[23]=1412.90215778 u[24]=1476.66111549 u[25]=1540.72049011 u[26]=1605.09292499 u[27]=1669.79107872 u[28]=1734.82762573 u[29]=1800.21525685 u[30]=1865.96667989 u[31]=1932.09462022 u[32]=1998.61182128 u[33]=2065.53104518 u[34]=2132.86507324 u[35]=2200.62670654 u[36]=2268.82876644 u[37]=2337.48409517 u[38]=2406.60555629 u[39]=2476.20603526 u[40]=2546.29843998 u[41]=2616.89570124 u[42]=2688.01077331 u[43]=2759.65663439 u[44]=2831.84628710 u[45]=2904.59275903 u[46]=2977.90910319 u[47]=3051.80839846 u[48]=3126.30375012 u[49]=3201.40829031 u[50]=3277.13517845 u[51]=3353.49760171 u[52]=3430.50877550 u[53]=3508.18194384 u[54]=3586.53037987 u[55]=3665.56738619 u[56]=3745.30629535 u[57]=3825.76047025 u[58]=3906.94330450 u[59]=3988.86822287 u[60]=4071.54868164 u[61]=4154.99816900 u[62]=4239.23020544 u[63]=4324.25834406 u[64]=4410.09617100 u[65]=4496.75730574 u[66]=4584.25540143 u[67]=4672.60414528 u[68]=4761.81725883 u[69]=4851.90849831 u[70]=4942.89165489 u[71]=5034.78055504 u[72]=5127.58906079 u[73]=5221.33107000 u[74]=5316.02051667 u[75]=5411.67137116 u[76]=5508.29764048 u[77]=5605.91336852 u[78]=5704.53263629 u[79]=5804.16956213 u[80]=5904.83830197 u[81]=6006.55304949 u[82]=6109.32803637 u[83]=6213.17753244 u[84]=6318.11584586 u[85]=6424.15732334 u[86]=6531.31635024 u[87]=6639.60735076 u[88]=6749.04478809 u[89]=6859.64316452 u[90]=6971.41702156 u[91]=7084.38094008 u[92]=7198.54954040 u[93]=7313.93748238 u[94]=7430.55946552 u[95]=7548.43022902 u[96]=7667.56455184 u[97]=7787.97725280 u[98]=7909.68319056 u[99]=8032.69726370 u[100]=8157.03441074 u[101]=8282.70961016 u[102]=8409.73788038 u[103]=8538.13427978 u[104]=8667.91390668 u[105]=8799.09189933 u[106]=8931.68343586 u[107]=9065.70373423 u[108]=9201.16805221 u[109]=9338.09168728 u[110]=9476.48997658 u[111]=9616.37829684 u[112]=9757.77206424 u[113]=9900.68673437 u[114]=10045.13780207 u[115]=10191.14080134 u[116]=10338.71130517 u[117]=10487.86492547 u[118]=10638.61731283 u[119]=10790.98415644 u[120]=10944.98118387 u[121]=11100.62416092 u[122]=11257.92889141 u[123]=11416.91121701 u[124]=11577.58701700 u[125]=11739.97220808 u[126]=11904.08274413 u[127]=12069.93461599 u[128]=12237.54385119 u[129]=12406.92651374 u[130]=12578.09870382 u[131]=12751.07655755 u[132]=12925.87624667 u[133]=13102.51397829 u[134]=13281.00599459 u[135]=13461.36857246 u[136]=13643.61802327 u[137]=13827.77069249 u[138]=14013.84295934 u[139]=14201.85123653 u[140]=14391.81196984 u[141]=14583.74163779 u[142]=14777.65675126 u[143]=14973.57385316 u[144]=15171.50951800 u[145]=15371.48035152 u[146]=15573.50299032 u[147]=15777.59410139 u[148]=15983.77038179 u[149]=16192.04855815 u[150]=16402.44538628 u[151]=16614.97765075 u[152]=16829.66216442 u[153]=17046.51576801 u[154]=17265.55532964 u[155]=17486.79774438 u[156]=17710.25993375 u[157]=17935.95884530 u[158]=18163.91145208 u[159]=18394.13475217 u[160]=18626.64576819 u[161]=18861.46154680 u[162]=19098.59915820 u[163]=19338.07569561 u[164]=19579.90827473 u[165]=19824.11403328 u[166]=20070.71013041 u[167]=20319.71374621 u[168]=20571.14208113 u[169]=20825.01235548 u[170]=21081.34180885 u[171]=21340.14769957 u[172]=21601.44730415 u[173]=21865.25791673 u[174]=22131.59684848 u[175]=22400.48142706 u[176]=22671.92899604 u[177]=22945.95691431 u[178]=23222.58255549 u[179]=23501.82330739 u[180]=23783.69657135 u[181]=24068.21976171 u[182]=24355.41030517 u[183]=24645.28564021 u[184]=24937.86321649 u[185]=25233.16049424 u[186]=25531.19494364 u[187]=25831.98404422 u[188]=26135.54528427 u[189]=26441.89616017 u[190]=26751.05417583 u[191]=27063.03684204 u[192]=27377.86167584 u[193]=27695.54619994 u[194]=28016.10794205 u[195]=28339.56443429 u[196]=28665.93321254 u[197]=28995.23181581 u[198]=29327.47778567 u[199]=29662.68866552 u[200]=30000.88200006
задание№4 #include <stdio.h> #include <math.h>
void main () { const double h = 0.2, xmax = 1.0, ymax = 1.0; const double eps = 0.01, Pi = 3.1415; const int kmax = 200; FILE * f; int n = (int) (xmax / h) + 1;
double ** u; u = new double * [n]; for (int i=0; i<n; i++) u[i] = new double [n]; fopen_s (&f, "station.dat", "w"); // Стороны AB и CD double y = 0.0; for (int j=0; j<n; j++) { u[0][j] =0.0; // Введите функцию на стороне AB u[n-1][j] = 50.0*y*(1-pow(y,2.0)); // Введите функцию на стороне CD y += h; } // Стороны AD и BC double x = 0.0; for (int i=0; i<n; i++) { u[i][0] = 50.0*x*(1-x); // Введите функцию на стороне AD u[i][n-1] =50.0*x*(1-x); // Введите функцию на стороне BC x += h; } // Начальные значения искомой функции внутри области for (int i=1; i<n-1; i++) for (int j=1; j<n-1; j++) u[i][j] = 0.0;
int k = 0; do { for (int i=1; i<n-1; i++) for (int j=1; j<n-1; j++) { u[i][j] = 0.25*(u[i-1][j]+u[i+1][j]+u[i][j+1]+u[i][j-1]); } k++; } while (k<kmax); fprintf (f, "Вычисленные значения\n"); for (int j=n-1; j>=0; j--) { for (int i=0; i<n; i++) { printf ("%2.4f\t", u[i][j]); fprintf (f, "%2.4f\t", u[i][j]); } printf ("\n"); fprintf (f, "\n"); }
for (int i=0; i<n; i++) delete [] u[i]; delete [] u; fclose (f); } Вычисленные значения 0.0000 8.0000 12.0000 12.0000 8.0000 0.0000 0.0000 5.3947 9.0948 10.8992 11.7296 14.4000 0.0000 4.4839 8.0853 10.7722 13.6192 19.2000 0.0000 4.4555 7.9904 10.4853 12.7748 16.8000 0.0000 5.3478 8.9355 10.4039 10.1947 9.6000 0.0000 8.0000 12.0000 12.0000 8.0000 0.0000
Задание№5
#include <stdio.h>
const int n = 3; const int n1 = 2;
int main () { double x[n] = { 0.0 }; double a[n][n] = { {15.5,-13.1,37.0}, {40.1,7.2,-11.9 },{11.5,54.8,23.7}}; //Введите коэффициенты при x, y, z double b[n] = {16.79,-53.23,141.57 }; // Введите столбец свободных членов
FILE * f; fopen_s (&f, "gauss.dat", "w");
// Печать исходной системы fprintf (f, "Исходная линейная система\n"); for (int l=0; l<n; l++) { for (int m=0; m<n; m++) { fprintf (f, "%2.3f * x[%d] ", a[l][m], m); if ((m<n-1) && (a[l][m+1]>=0.0)) fprintf (f, "+"); if (m==(n-1)) fprintf (f, " = "); } fprintf (f, "%2.3f\n", b[l]); }
// Прямой ход метода Гаусса double w[n1]; for (int i=1; i<n; i++) { int l = i - 1; for (int j=0; j<n-i; j++) w[j] = -a[j+i][l] / a[l][l]; for (int k=i; k<n; k++) { for (int j=0; j<n; j++) a[k][j] += a[l][j]*w[k-i]; b[k] += b[l]*w[k-i]; } }
// Печать системы, приведенной к треугольному виду fprintf (f, "\nСистема, приведенная к треугольному виду\n"); for (int l=0; l<n; l++) { for (int m=0; m<n; m++) { fprintf (f, "%2.3f * x[%d] ", a[l][m], m); if ((m<n-1) && (a[l][m+1]>=0.0)) fprintf (f, "+"); if (m==(n-1)) fprintf (f, " = "); } fprintf (f, "%2.3f\n", b[l]); }
// Обратный ход метода Гаусса double sum; for (int j=n-1; j>=0; j--) { sum = 0.0; for (int i=j; i<n-1; i++) sum += a[j][i+1]*x[i+1]; x[j] = (b[j] - sum) /a[j][j]; }
// Печать результатов printf ("Solution - \n"); fprintf (f, "\nВычисленные решения системы\n"); for (int i=0; i<n; i++) { printf ("x[%d] = %2.4f ", i, x[i]); fprintf (f, "x[%d] = %2.4f ", i, x[i]); } printf ("\n"); fprintf (f, "\n");
fclose (f); return 0; }
Исходная линейная система 15.500 * x[0] -13.100 * x[1] +37.000 * x[2] = 16.790 40.100 * x[0] +7.200 * x[1] -11.900 * x[2] = -53.230 11.500 * x[0] +54.800 * x[1] +23.700 * x[2] = 141.570
Система, приведенная к треугольному виду 15.500 * x[0] -13.100 * x[1] +37.000 * x[2] = 16.790 0.000 * x[0] +41.091 * x[1] -107.623 * x[2] = -96.667 0.000 * x[0] +0.000 * x[1] +165.233 * x[2] = 280.896
Вычисленные решения системы x[0] = -1.2000 x[1] = 2.1000 x[2] = 1.7000
|
