الإنجليزية: Implicit Euler Method
الفرنسية: Méthode d'Euler Implicite
المبدأ
حل رقميا y'(x) = f(x, y(x)) x [a, b] y(a) = y0
نجزئ المجال الاندماج[a, b] إلىN مجال ذي طول يساوي h = ،
محدود بنقط التجزيء xn = a + nh , n = 0,... N .
لكل أفصول xn، نحسب القيمة المقربة yn للقيمة الصحيحة y(xn) لدالة y .
لهدا نعتبر مستقيم ميل ( f(xn+1, yn+1 المار من نقطة (xn, yn) (موازي لمماس المنحنى التام المار من النقطة (xn+1, yn+1))
القيمة yn+1 هي أرتوب نقطة هدا المستقيم ذي أفصول xn+1.
الخوارزم
(إنها معادلة من yn+1 الواجب حلها باستعمال الطريقة الرقمية)
مثال
نعتبر مشكلة كوشي:
لتكن علما أن y(0) = 1
الحل المضبوط هو :
احسب رقميا قيمy(x) لـ x = 0, 0.1, 0.2, ..., 1 بأخد لـ N القوى المتتالية لـ 10.
الخطوتين الأولتين لـ N = 10
h = 1/10
x0 = 0 , y0 = 1
x1 = 1/10
y1 = y0 + h f(x0+h, y1)
= y0 - h y12
= 1 - (1/10)y12
y1 هو حل المعادلة
y12 + 10y1 - 10 = 0
ومنه
y1= -5 + ?35 |
x2 = 2/10
y2 = y1 + h f(x1+h, y2) |
y2 هو حل المعادلة
y22 + 10y2 + 50 - 10√35 = 0 |
ومنه
y2= -5 + |
البرمجة
في كل خطوة، المعادلة المحصل عليها نحلها باستعمال طريقة ستفنسن
#include <math.h> #include <conio.h> #define ITERMAX 11 double x,y,h;
double ed_f(double x,double y) { return(-y*y); }
double eq_f(double yy); double eq_g(double x); void eq_steffensen(double x,double eps,double t[ITERMAX]);
double ed_solution(double x) { return(1.0/(1.0+x)); }
long ed_euler_implicite(double a,double b,double y0,long n,double res[101][2]); void ed_affichage(double a,double b,long n,long nmax,double res[101][2],int nbpts);
main() { int nbpts; long n,nmax; double y0,a,b; double res[101][2]; clrscr(); nbpts=10; a=0;b=1;y0=1; printf("Méthode d'Euler implicite\n"); for(n=10;n<=10000;n*=10) { nmax=ed_euler_implicite(a,b,y0,n,res); ed_affichage(a,b,n,nmax,res,nbpts); printf(" taux = %.3e\n", fabs(res[nmax][1]-ed_solution(res[nmax][0]))*n); } }
نتيجة البرنامج
-
في كل خطوة البرنامج يظهر:
-
قيمة xn
-
قيمة yn
-
الارتياب (الفرق بين yn المحسوبة والقيمة المضبوطةy(xn) )
ويظهر أيضا، لكل قيمة h، نسبة التقارب، العلاقة بين الارتياب وh في نقطة اعتباطية من المجال (نختار x = 1)
طريقة ؤلير الضمني
h = 0.1000 xn yn الارتياب 0.0000 1.0000000000000000e+00 0.00e+00 0.1000 9.1607978309961602e-01 6.99e-03 0.2000 8.4472393111908750e-01 1.14e-02 0.3000 7.8335882607943275e-01 1.41e-02 0.4000 7.3006005734619916e-01 1.58e-02 0.5000 6.8336173170967462e-01 1.67e-02 0.6000 6.4212879302632953e-01 1.71e-02 0.7000 6.0546946564364901e-01 1.72e-02 0.8000 5.7267392339050194e-01 1.71e-02 0.9000 5.4317050377354170e-01 1.69e-02 1.0000 5.1649390806655537e-01 1.65e-02 نسبة = 1.65e-01 h = 0.0100 xn yn الارتياب 0.0000 1.0000000000000000e+00 0.00e+00 0.1000 9.0986851721306483e-01 7.78e-04 0.2000 8.3458523410850904e-01 1.25e-03 0.3000 7.7076777107117478e-01 1.54e-03 0.4000 7.1598715042494554e-01 1.70e-03 0.5000 6.6845432596457455e-01 1.79e-03 0.6000 6.2682266589573776e-01 1.82e-03 0.7000 5.9005928539682384e-01 1.82e-03 0.8000 5.5735878334768885e-01 1.80e-03 0.9000 5.2808393731517855e-01 1.77e-03 1.0000 5.0172401987026161e-01 1.72e-03 نسبة = 1.72e-01 h = 0.0010 xn yn الارتياب 0.0000 1.0000000000000000e+00 0.00e+00 0.1000 9.0916957564373957e-01 7.87e-05 0.2000 8.3345980157761124e-01 1.26e-04 0.3000 7.6938585818388849e-01 1.55e-04 0.4000 7.1445722971314851e-01 1.72e-04 0.5000 6.6684672802478551e-01 1.80e-04 0.6000 6.2518346123111768e-01 1.83e-04 0.7000 5.8841878067665121e-01 1.83e-04 0.8000 5.5573686114564713e-01 1.81e-04 0.9000 5.2649348928419393e-01 1.78e-04 1.0000 5.0017319776946068e-01 1.73e-04 نسبة = 1.73e-01 h = 0.0001 xn yn الارتياب 0.0000 1.0000000000000000e+00 0.00e+00 0.1000 9.0909878494194318e-01 7.88e-06 0.2000 8.3334599311131308e-01 1.27e-05 0.3000 7.6924629218001706e-01 1.55e-05 0.4000 7.1430287969423956e-01 1.72e-05 0.5000 6.6668468588327445e-01 1.80e-05 0.6000 6.2501835817630924e-01 1.84e-05 0.7000 5.8825365373754268e-01 1.84e-05 0.8000 5.5557369601804485e-01 1.81e-05 0.9000 5.2633356836857814e-01 1.78e-05 1.0000 5.0001732778870567e-01 1.73e-05 نسبة = 1.73e-01
الدالة الرئيسية ed_euler_implicite
معايير هذه الدالة هي كالتالي:
-
a و b : حدي المجال [a, b]
-
y0 : الشرط البدئي y0 = y(a) = y(x0)i
-
n : عدد المجالات الصغيرة (ذات الطول h)
-
res : جدول لحفظ القيم (xn, yn)
تستدعي الدالة الرئيسية الدالة التالية التي يجب تعريفها في البرنامج الرئيسي أيضا:
double ed_f(double x,double y)
تمثل هذه الدالة الدالة الرياضية f .
وتستدعي أيضا دالة ستفنسن المعرفة كالتالي:
void eq_steffensen(double x,double eps,double t[ITERMAX])
هذه الدالة بدورها تستدعي دالتين أخريين:
يجب أن تعرف جميع هذه الدوال في البرنامج الرئيسي.
وأن قيمة الثابتة الصحيحة ITERMAX تساوي العدد القصوي لمرات التكرار الخاصة بطريقة ستفنسن (+1).
المتغيرات الثلاث x و y و h يجب أن تعرف كمتغيرات عامة في البرنامج.
يتم حساب القيم التقريبية لـ y(x)i من خلال طريقة أولير الضمنية. وترجع بعدد القيم (xn, yn) المحفوظة (أقصى حد هو 100).
long ed_euler_implicite(double a,double b,double y0,long n,double res[101][2])
{
double tol,xx,exx;
long i,j,k;
double t[ITERMAX];
y=y0;
h=(b-a)/n;
tol=h/2;
x=a;
res[0][0]=a;
res[0][1]=y; j=1;
for(i=1;i<=n;i++)
{
eq_steffensen(y+h*ed_f(x,y),1e-14,t);
k=ITERMAX-1;
while(t[k]==0) k--;
y=t[k];
x += h;
xx=x-a;
exx=floor(100*xx+0.5)/100;
if(fabs(exx-xx)<tol)
{
res[j][0]=x;
res[j][1]=y;
j++;
}
}
return(j-1);
}
طريقة ستفنسن eq_steffensen
معايير هذه الدالة هي كالتالي:
- x : يمثل x0 للانطلاق
- eps :الدقة المطلوبة
- t :جدول ذي البعد ITERMAX.
تستدعي هذه الدالة الدالتين:
double eq_f(double x)
double eq_g(double x)
ترجع هذه الدالة في الجدول t التكرارات المتتابعة حسب طريقة ستفنسن.
قيمة الثابتة الصحيحة ITERMAX تساوي العدد القصوي لمرات التكرار الخاصة بطريقة ستفنسن (+1).
void eq_steffensen(double x,double eps,double t[ITERMAX])
{
double y,err,den1,den2,den3;
int n,stop;
for(n=0;n<ITERMAX;n++) t[n]=0;
n=0;
err=2*eps;
stop=0;
while ((fabs(err) > eps) && (stop==0) && (n < ITERMAX-1))
{
n++;
den1=eq_g(eq_g(x))-eq_g(x);
den2=eq_g(x)-x;
if((den1!=0) && (den2!=0))
{
den3=1.0/den1-1.0/den2;
if(den3!=0)
{
y=eq_g(x)+1.0/den3;
err=y-x;
x=y;
t[n]=x;
}
else
stop=1;
}
else stop=1;
}
}
دالة الإظهار ed_affichage
معايير دالة الإظهار هي كالتالي:
-
a و b : حدي المجال [a, b]
n : عدد المجالات الصغيرة (ذات الطول h)
nmax : العدد الأقصى للقيم المحتفظ بها
res : جدول يحفظ القيم (xn, yn)
nbpts : عدد القيم المسموح لها بالإظهار على الشاشة
وهي بدورها تستدعي الدالة التالية:
double ed_solution(double x)
التي يجب أن تعرف في البرنامج الرئيسي وتعطى الحل الحقيقي لـ y.
ستقوم دالة الإظهار بإظهار جزئي للقيم (nbpts قيمة كأقصى حد) المحفوظة في الجدول res. وتظهر أيضا نسبة الخطأ، أي الفرق بينالقيمة الحقيقية y(xn)i للحل y(x)i والقيمة المقربة لـ yn.
void ed_affichage(double a,double b,long n,long nmax,double res[101][2],int nbpts)
{
double h,x;
long i;
h=(b-a)/n;
printf("h = %.4f\n",h);
printf(" xn yn Erreur\n");
if(nmax<nbpts)
for(i=0;i<=nmax;i++)
printf("%8.4f %.17e %.3e\n", res[i][0],res[i][1],fabs(ed_solution(res[i][0])-res[i][1]));
else
{
for(i=0;i<=nmax;i++)
{
x=floor(res[i][0]*nbpts+0.5)/nbpts;
if(fabs(x-res[i][0])<0.005)
printf("%8.4f %.17e %.3e\n", x,res[i][1],fabs(ed_solution(res[i][0])-res[i][1]));
}
}
}