أجديك عبد الرحيم http://www.isnaha.com Wed, 18 Jul 2018 12:41:08 +0000 Joomla! - Open Source Content Management ar-aa التقريب الرياضي http://www.isnaha.com/isnaha_new/برمجها/item/962-التقريب-الرياضي http://www.isnaha.com/isnaha_new/برمجها/item/962-التقريب-الرياضي

التقريب الرياضي

 approxi
 

الإشكالية
 
نعتبر : دالة متصلة أو غير متصلة. 
لنحدد دالة أخرى g :  ذات صيغة معطاة بحيث تقترب من الدالة f من أحد الاتجاهات. 

الخوارزميات 
 
نفرد لحل هذه الإشكالية خوارزميتين (إضغط على الرابط):

الخوارزم 1: التقريب الجذري

الخوارزم 2: التقريب في اتجاه صغار المربعات


]]>
aj.abderahim@gmail.com (أجديك عبد الرحيم) التقريب Mon, 24 Jun 2013 07:23:17 +0000
طريقة الفوارق المنتهية http://www.isnaha.com/isnaha_new/برمجها/item/954-طريقة-الفوارق-المنتهية http://www.isnaha.com/isnaha_new/برمجها/item/954-طريقة-الفوارق-المنتهية

طريقة الفوارق المنتهية

differences-finies 
 

مصطلحات
العربية: طريقة الفوارق المنتهية
الإنجليزية: Finite Difference Method
الفرنسية: Méthode des Différences Finies

المبدأ

 حل رقميا   y'(x) = f(xy(x))    x  [ab]   y(a) = y0

نجزئ المجال الاندماج [a, b] إلى N مجال ذي طول يساوي h = ،
محدودة بنقط التجزيء xn = a + nh , n = 0,... N .

لكل أفصول xn، نحسب القيمة المقربة yn للقيمة الصحيحة y(xn) لدالة y .

خاصية طايلور Taylor تمكننا من كتابة :

بتعويض في المعادلة التفاضلية بتقريب المحصل عليه بإهمال لعبارات من h2.


الخوارزم

تحل النظمة الخطية ذات N-1 معادلة وN-1  مجهول .

       n = 1,... N-1

(yN وy0 معطات)


مثال

أجرى عملية التكامل :

 x [0, 1]  بحيث      علما أن y(0)=1 و  y(1)=1/2  

الحل المضبوط هو : 

احسب رقميا قيمy(x)   لـ x = 0, 0.1, 0.2, ..., 1  بأخد  لـ N القوى المتتالية لـ 10.

النظمة لـ N = 10

h = 1/10
x0 = 0 , y0 = 1 ; x10 = 1 , y10 = 1/2

المعادلة الأولى

(y0-2y1+y2)/(10-2) = f(x1) = f(1/10) = 2/(1+1/10)3
 y0 - 2y1 + y2 = 20/1331       ومنه

وبما أن y0 = 1

-2y1 + y2 = -1311/1331

المعادلة الثانية

(y1-2y2+y3)/(10-2) = f(x2) = f(2/10) = 2/(1+2/10)3

 y1 - 2y2 + y3 = 5/432     ومنه

وبالمثل للمعادلات الأخرى

y2 - 2y3 + y4 = 20/2197
y3 - 2y4 + y5 = 5/686
y4 - 2y5 + y6 = 4/675
y5 - 2y6 + y7 = 5/1024
y6 - 2y7 + y8 = 20/4913
y7 - 2y8 + y9 = 5/1458
y8 - 2y9 = -6819/13718


البرمجة

#include <math.h>
#include <conio.h>
#define NMAX 1002
#define NDIAG 2

void al_prod_matmd_vect(double a[NMAX][NDIAG],long n,double x[NMAX],double y[NMAX]); void sl_grad_conj(double a[NMAX][NDIAG],double b[NMAX],double x1[NMAX],long n);
double ed_f(double x) { return(2.0/pow(1.0+x,3)); }
double ed_solution(double x) { return(1.0/(1.0+x)); }
long ed_differences_finies(double a,double b,double y0,double yn,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 i,n,nmax; double a,b,y0,yn; double res[101][2]; clrscr(); nbpts=10; a=0;b=1;y0=1;yn=0.5; printf("Méthode des différences finies\n"); for(n=10;n<=1000;n*=10) { nmax=ed_differences_finies(a,b,y0,yn,n,res); ed_affichage(a,b,n,nmax,res,nbpts); for(i=1;i<=nmax;i++) if(fabs(res[i][0]-0.5)*n<0.5) printf(" taux = %.3e\n", fabs(res[i][1]-ed_solution(res[i][0]))*n*n); } }

نتيجة البرنامج

  • في كل خطوة البرنامج يظهر:
    قيمة xn
    قيمة yn
    الارتياب (الفرق بين yالمحسوبة والقيمة المضبوطةy(xn) )

ويظهر أيضا، لكل قيمة h، نسبة التقارب، العلاقة بين الارتياب وh في نقطة اعتباطية من المجال (نختار x = 0.5).

h =  0.1000
     xn        yn              الارتياب
  0.0000  1.0000000000000000e+00  0.00e+00
  0.1000  9.0935677412321680e-01  2.66e-04
  0.2000  8.3373984426446524e-01  4.07e-04
  0.3000  7.6969698847978774e-01  4.66e-04
  0.4000  7.1475745540790048e-01  4.72e-04
  0.5000  6.6710655207362257e-01  4.40e-04
  0.6000  6.2538157466527067e-01  3.82e-04
  0.7000  5.8853940975691876e-01  3.04e-04
  0.8000  5.5576807733381006e-01  2.13e-04
  0.9000  5.2642610019190861e-01  1.10e-04
  1.0000  5.0000000000000000e-01  2.78e-17
      نسبة = 4.40e-02
h =  0.0100
     xn        yn              الارتياب
  0.0000  1.0000000000000000e+00  0.00e+00
  0.1000  9.0909359522873889e-01  2.69e-06
  0.2000  8.3333743786729886e-01  4.10e-06
  0.3000  7.6923547435519712e-01  4.71e-06
  0.4000  7.1429047333442908e-01  4.76e-06
  0.5000  6.6667110301117549e-01  4.44e-06
  0.6000  6.2500384733290926e-01  3.85e-06
  0.7000  5.8823835983916306e-01  3.07e-06
  0.8000  5.5555769758743101e-01  2.14e-06
  0.9000  5.2631690115516461e-01  1.11e-06
  1.0000  5.0000000000000000e-01  1.67e-16
      نسبة = 4.44e-02
h =  0.0010
     xn       yn              الارتياب
  0.0000  1.0000000000000000e+00  0.00e+00
  0.1000  9.0909093595506119e-01  2.69e-08
  0.2000  8.3333337438264432e-01  4.10e-08
  0.3000  7.6923081628636714e-01  4.71e-08
  0.4000  7.1428576188051895e-01  4.76e-08
  0.5000  6.6666671103405495e-01  4.44e-08
  0.6000  6.2500003847669205e-01  3.85e-08
  0.7000  5.8823532477751395e-01  3.07e-08
  0.8000  5.5555557697772795e-01  2.14e-08
  0.9000  5.2631580059145200e-01  1.11e-08
  1.0000  5.0000000000000000e-01  1.67e-16
     نسبة = 4.44e-02

الدالة الرئيسية ed_euler

 

معايير الدالة al_prod_mat_vect هي كالتالي:

a : المصفوفة
n : درجة المصفوفة
x و y : متجهتين

ترجع هذه الدالة في المتجهة y ضرب المصفوفة a في المتجهة x.

 الدالة al_prod_matmd_vect

 المصفوفة متماثلة قطريا، ومحفوظة في جدول ثنائي البعد بعده NMAX*NDIAG.

 معايير الدالة هي كما يلي :

 a : مصفوفة متماثلة قطريا ومكدسة

 n : درجة المصفوفة

 x و y : متجهتين

 ترجع هذه الدالة في المتجهة y ضرب المصفوفة a في المتجهة x.

الدالة sl_grad_conj
معاييرها هي كالتالي :
a : مصفوفة النظام
b : جدول أحادي البعد (يمثل العنصر الثانس)
x1 :جدول أحادي البعد (متجهة الإنطلاق) 
n : درجة المصفوفة
t : مصفوفة

تستدعي هذه الدالة بدورها الدالة: 

void al_prod_mat_vect(double a[NMAX][NMAX],int n,double x[NMAX],double y[NMAX])

التي يجب أن تعرف في البرنامج الرئيسي.

ترجع الدالة الرئيسية في الجدول t التكرارات المتتالية لطريقة الفوارق المنتهية ; تملئ الجدول انطلاقا من الخانات الأكثر قبولا للقيم الجديدة.

الثابتة الصحيحة NMAX تساوي البعد القصوي للمصفوفة (+1).

void sl_grad_conj(double a[NMAX][NMAX],double b[NMAX],double x1[NMAX], int n,double t[NMAX][NMAX])
{
int i,j,k;
double alfa0,alfa1,s,mu,lambda;
double r[NMAX],u[NMAX],v[NMAX];

for(i=1;i<=n;i++)
for(j=1;j<=NMAX-1;j++) t[i][j]=0.0;

al_prod_mat_vect(a,n,x1,r);
alfa0=0;
for(i=1;i<=n;i++)
{
r[i] -= b[i];
u[i]=-r[i];
alfa0 += r[i]*r[i];
}
for(k=0;k<=n;k++)
{
al_prod_mat_vect(a,n,u,v);
s=0;
for(i=1;i<=n;i++) s += u[i]*v[i];
mu=alfa0/s;
for(i=1;i<=n;i++)
{
x1[i] += mu*u[i];
for(j=1;j<=NMAX-2;j++) t[i][j]=t[i][j+1];
t[i][NMAX-1]=x1[i];
}
al_prod_mat_vect(a,n,x1,r);
alfa1=0;
for(i=1;i<=n;i++)
{
r[i] -= b[i];
alfa1 += r[i]*r[i];
}
if (alfa1==0) break;
lambda=alfa1/alfa0;
alfa0=alfa1;
for(i=1;i<=n;i++)
u[i]=-r[i]+lambda*u[i];
}
}

 دالة الفوارق المنتهية ed_differences_finies


معايير هذه الدالة هي كالتالي:

 a و b : حدي المجال [ab]
y0 : الشرط البدئي y0 = y(a) = y(x0)i
yn : الشرط البدئي yN = y(b) = y(xN)i
n : عدد المجالات الصغيرة (ذات الطول h)
res : جدول لحفظ القيم (xnyn)

تستدعي الدالة الرئيسية الدالة التالية التي يجب تعريفها في البرنامج الرئيسي أيضا:

 double ed_f(double x,double y)

تمثل هذه الدالة الدالة الرياضية  f .
وتستدعي أيضا الدالة:

void sl_grad_conj(double a[NMAX][NDIAG],double b[NMAX], double x1[NMAX], long n)

يتم حساب القيم التقريبية لـ y(x)i من خلال طريقة الفوارق المنتهية. وترجع بعدد القيم (xnyn) المحفوظة (أقصى حد هو 100). سيتم حفظ هذه القيم في الجدول res :حيث سيحتوي العمود 0 على قيم xn بينما سيمثل العمود 1 قيم yn المقابلة.

 قيمة الثابتة NDIAG تساوي 2.

long ed_differences_finies(double a,double b,double y0,double yn,long n,double res[101][2]) 
{
double x,h,tol,xx,exx;
double m[NMAX][NDIAG],g[NMAX],y[NMAX];
long i,j;
h=(b-a)/n;
for(i=1;i<n;i++)
{
m[i][0]=2.0/(h*h);
if(i<n-1)m[i][1]=-1.0/(h*h);
}
x=a;
for(i=1;i<n;i++)
{
x+=h;
g[i]=-ed_f(x);
y[i]=1.0;
}
g[1]+=y0/(h*h);
g[n-1]+=yn/(h*h);
sl_grad_conj(m,g,y,n-1);
y[n]=yn;
tol=h/2;
x=a;
res[0][0]=a;
res[0][1]=y0;
j=1;
for(i=1;i<=n;i++)
{
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[i];
j++;
}
}
return(j-1);
}

 دالة الإظهار ed_affichage

معايير دالة الإظهار هي كالتالي:

  • a و b : حدي المجال [ab]
    n : عدد المجالات الصغيرة (ذات الطول h)
    nmax : العدد الأقصى للقيم المحتفظ بها
    res : جدول يحفظ القيم (xnyn)
    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]));
}

}
}
]]>
aj.abderahim@gmail.com (أجديك عبد الرحيم) التكامل الرقمي للمعادلات التفاضلية العادية Tue, 18 Jun 2013 13:33:49 +0000
طريقة آدم باشفورت مولتون http://www.isnaha.com/isnaha_new/برمجها/item/953-طريقة-آدم-باشفورت-مولتن http://www.isnaha.com/isnaha_new/برمجها/item/953-طريقة-آدم-باشفورت-مولتن

طريقة آدم باشفورت مولتون

 adam b m
 

مصطلحات
العربية: طريقة آدم باشفورت مولتون
الإنجليزية: Adams Bashforth Moulton Method
الفرنسية: Méthode d'Adams Bashforth Moulton

المبدأ

 حل رقميا   y'(x) = f(xy(x))    x  [ab]   y(a) = y0

 نجزئ المجال الاندماج [a, b] إلى N مجال ذي طول يساوي h = ،
محدودة بنقط التجزيء xn = a + nh , n = 0,... N .

لكل أفصول xn، نحسب القيمة المقربةyn  للقيمة الصحيحة y(xn) لدالة y .

لهذا، نستعمل الخاصية التالية:

نقرب  على  بحدودية الاستكمال في النقط:

  • xn-3, xn-2, xn-1, xn (صيغة التنبؤ)
  • xn-2, xn-1, xn, xn+1 (صيغة التصحيح)
بتعويض  y(xn)بـ yn لكل n.

الخوارزم

تهيئ:

حسابy1 و y2 و y3 باستعمال طريقة ذات خطوة واحدة

الكرة:

صيغة التنبؤ:

صيغة التصحيح:


مثال

نعتبر مشكلة كوشي:

لتكن   علما أن y(0) = 1

الحل المضبوط هو : 

احسب رقميا قيمy(x)   لـ x = 0, 0.1, 0.2, ..., 1  بأخد  لـ N القوى المتتالية لـ 10.

 3 خطوات الأولى (حساب y1, y2, y3) سيتم باستعمال طريقة رونج-كيتا.

الخطوتين الأولتين لـ N = 10

h = 1/10
x0 = 0 , y0 = 1
x1 = 1/10 , x2 = 2/10 , x3 = 3/10 ( y1, y2, y3 و منه )

x4 = 4/10
y4* = y3 + (1/240) [55f(x3, y3) - 59f(x2, y2) + 37f(x1, y1) - 9f(x0, y0)]
y4 = y3 + (1/240) [9f(x4, y4*) + 19f(x3, y3) - 5f(x2, y2) + f(x1, y1)]


البرمجة

#include <math.h>
#include <conio.h>
double ed_f(double x,double y) { return(-y*y); }
double ed_solution(double x) { return(1.0/(1.0+x)); }
long ed_adams_bashforth_moulton(double a,double b,double y,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 a,b,y0; double res[101][2]; clrscr(); nbpts=10; a=0;b=1;y0=1; printf("Méthode d'Adams-Bashforth-Moulton\n"); for(n=10;n<=1000;n*=10) { nmax=ed_adams_bashforth_moulton(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]))*pow(n,4)); } }

نتيجة البرنامج

  • في كل خطوة البرنامج يظهر:
    قيمة xn
    قيمة yn
    الارتياب (الفرق بين yالمحسوبة والقيمة المضبوطةy(xn) )

ويظهر أيضا، لكل قيمة h، نسبة التقارب، العلاقة بين الارتياب وh4 في نقطة اعتباطية من المجال (نختار x = 1) 

h =  0.1000
     xn             yn            الارتياب
  0.0000  1.0000000000000000e+00  0.00e+00
  0.1000  9.0909118633221964e-01  2.77e-07
  0.2000  8.3333372884307211e-01  3.96e-07
  0.3000  7.6923120575328652e-01  4.37e-07
  0.4000  7.1427030768034339e-01  1.54e-05
  0.5000  6.6664403261030503e-01  2.26e-05
  0.6000  6.2497434530089857e-01  2.57e-05
  0.7000  5.8820877835436980e-01  2.65e-05
  0.8000  5.5552937890959886e-01  2.62e-05
  0.9000  5.2629056824798792e-01  2.52e-05
  1.0000  4.9997602892674858e-01  2.40e-05
         نسبة = 2.40e-01
h =  0.0100
     xn             yn            الارتياب
  0.0000  1.0000000000000000e+00  0.00e+00
  0.1000  9.0909090744887688e-01  1.64e-09
  0.2000  8.3333333050683434e-01  2.83e-09
  0.3000  7.6923076594958040e-01  3.28e-09
  0.4000  7.1428571090882076e-01  3.38e-09
  0.5000  6.6666666336955593e-01  3.30e-09
  0.6000  6.2499999686463836e-01  3.14e-09
  0.7000  5.8823529117751572e-01  2.94e-09
  0.8000  5.5555555281893010e-01  2.74e-09
  0.9000  5.2631578693596093e-01  2.54e-09
  1.0000  4.9999999765032532e-01  2.35e-09
       نسبة = 2.35e-01
h =  0.0010
     xn             yn            الارتياب
  0.0000  1.0000000000000000e+00  0.00e+00
  0.1000  9.0909090909069668e-01  2.12e-13
  0.2000  8.3333333333302617e-01  3.07e-13
  0.3000  7.6923076923042988e-01  3.39e-13
  0.4000  7.1428571428537246e-01  3.42e-13
  0.5000  6.6666666666633689e-01  3.30e-13
  0.6000  6.2499999999968892e-01  3.11e-13
  0.7000  5.8823529411735709e-01  2.90e-13
  0.8000  5.5555555555528657e-01  2.69e-13
  0.9000  5.2631578947343471e-01  2.49e-13
  1.0000  4.9999999999976941e-01  2.30e-13
   نسبة = 2.30e-01

الدالة الرئيسية ed_adams_bashforth_moulton

معايير هذه الدالة هي كالتالي:

  • a و b : حدي المجال [ab]

  • y : الشرط البدئي y0 = y(a) = y(x0)i

  • n : عدد المجالات الصغيرة (ذات الطول h)

  • res : جدول لحفظ القيم (xnyn)

تستدعي الدالة الرئيسية الدالة التالية التي يجب تعريفها في البرنامج الرئيسي أيضا:

 double ed_f(double x,double y)

 تمثل هذه الدالة الدالة الرياضية  f .

يتم حساب القيم التقريبية لـ y(x)i من خلال طريقة آدم باشفورت مولتون. وترجع بعدد القيم (xnyn) المحفوظة (أقصى حد هو 100). سيتم حفظ هذه القيم في الجدول res: حيث سيحتوي العمود 0 على قيم xn بينما سيمثل العمود 1 قيم yn المقابلة.

ملاحظة: الخطوات الثلاث الأولى تتم بمساعدة طريقة رونج كيتا بدرجة 4.

long ed_adams_bashforth_moulton(double a,double b,double y,long n,double res[101][2]) { double x,h,k1,k2,k3,k4,tol,xx,yy,exx; double f[4]; long i,j,k; 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++) { if(i==1||i==2||i==3) { k1=ed_f(x,y); k2=ed_f(x+h/2,y+h*k1/2); k3=ed_f(x+h/2,y+h*k2/2); k4=ed_f(x+h,y+h*k3); f[i]=k1; y += h*(k1+2*k2+2*k3+k4)/6; } else { for(k=0;k<=2;k++) f[k]=f[k+1]; f[3]=ed_f(x,y); yy = y+h*(55*f[3]-59*f[2]+37*f[1]-9*f[0])/24; y += h*(9*ed_f(x+h,yy)+19*f[3]-5*f[2]+f[1])/24; } 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); } 

 دالة الإظهار ed_affichage

معايير دالة الإظهار هي كالتالي:

  • a و b : حدي المجال [ab]
    n : عدد المجالات الصغيرة (ذات الطول h)
    nmax : العدد الأقصى للقيم المحتفظ بها
    res : جدول يحفظ القيم (xnyn)
    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]));
}

}
}
]]>
aj.abderahim@gmail.com (أجديك عبد الرحيم) التكامل الرقمي للمعادلات التفاضلية العادية Tue, 18 Jun 2013 13:32:44 +0000
طريقة أولير الضمنية http://www.isnaha.com/isnaha_new/برمجها/item/952-طريقة-أولير-الضمنية http://www.isnaha.com/isnaha_new/برمجها/item/952-طريقة-أولير-الضمنية
طريقة أولير الضمنية
 Euler-implicite
 

مصطلحات
العربية: طريقة أولير الضمنية
الإنجليزية: Implicit Euler Method
الفرنسية: Méthode d'Euler Implicite

المبدأ

حل رقميا   y'(x) = f(xy(x))    x  [ab]   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)
     = y1 - h y22
     = -5 + √35 - (1/10)y22

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

  • الارتياب (الفرق بين yالمحسوبة والقيمة المضبوطة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 : حدي المجال [ab]

  • y0 : الشرط البدئي y0 = y(a) = y(x0)i

  • n : عدد المجالات الصغيرة (ذات الطول h)

  • res : جدول لحفظ القيم (xnyn)

تستدعي الدالة الرئيسية الدالة التالية التي يجب تعريفها في البرنامج الرئيسي أيضا:

 double ed_f(double x,double y)

 تمثل هذه الدالة الدالة الرياضية  f .

وتستدعي أيضا دالة ستفنسن المعرفة كالتالي:

void eq_steffensen(double x,double eps,double t[ITERMAX]) 

هذه الدالة بدورها تستدعي دالتين أخريين:

double eq_f(double yy)
double eq_g(double x)

يجب أن تعرف جميع هذه الدوال في البرنامج الرئيسي.
وأن قيمة الثابتة الصحيحة ITERMAX تساوي العدد القصوي لمرات التكرار الخاصة بطريقة ستفنسن (+1).

المتغيرات الثلاث x  و y و  h يجب أن تعرف كمتغيرات عامة في البرنامج.

يتم حساب القيم التقريبية لـ y(x)i من خلال طريقة أولير الضمنية. وترجع بعدد القيم (xnyn) المحفوظة (أقصى حد هو 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 : حدي المجال [ab]
    n : عدد المجالات الصغيرة (ذات الطول h)
    nmax : العدد الأقصى للقيم المحتفظ بها
    res : جدول يحفظ القيم (xnyn)
    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]));
}

}
}
]]>
aj.abderahim@gmail.com (أجديك عبد الرحيم) التكامل الرقمي للمعادلات التفاضلية العادية Tue, 18 Jun 2013 13:31:48 +0000
طريقة نيستروم http://www.isnaha.com/isnaha_new/برمجها/item/951-طريقة-نيستروم http://www.isnaha.com/isnaha_new/برمجها/item/951-طريقة-نيستروم
طريقة نيستروم
 nystrom
 

مصطلحات

العربية: طريقة نيستروم
الإنجليزية: Nyström Method
الفرنسية: Méthode de Nyström
 


المبدأ

 حل رقميا   y'(x) = f(xy(x))    x  [ab]   y(a) = y0

نجزئ المجال الاندماج [a, b] إلى N مجال ذي طول يساوي h = ،
محدودة بنقط التجزيء xn = a + nh , n = 0,... N .

لكل أفصول xn، نحسب القيمة المقربة yn للقيمة الصحيحة y(xn) لدالة y.

لهدا نعتبر مستقيم ميل ( f(xn, yn المار من نقطة(xn-1, yn-1) (موازي  لمماس المنحنى التام  المار من النقطة (xn, yn)) القيمةyn+1 هي أرتوب نقطة هدا المستقيم ذي أفصولxn+1.


مثال

نعتبر مشكلة كوشي: 
مع y(0) = 1

الحل المضبوط هو :

احسب رقميا قيمy(x)   لـ x = 0, 0.1, 0.2, ..., 1  بأخد  لـ N القوى المتتالية لـ 10.

الخطوة الأولى(حساب y1) سيتم باستعمال طريقة ؤلير المتطورة.

الخطوتين الأولتين لـ N = 10

h = 1/10
x0 = 0 , y0 = 1

x1 = 1/10

الخطوة الأولى طريقة ؤلير المتطورة

y1 = y0 + h f(x0+(h/2), y0+(h/2) f(x0, y0))
     = y0 - h( y0 -(h/2)y02 )2
     = 1 - (1/10)(1-(1/20)12)2
     = 3639/4000

x2 = 2/10
y2 = y0 + 2h f(x1, y1)
     = y0 - 2h y12
     = 1 - 2(1/10)(3639/4000)2
     = 66757679/80000000


البرمجة

#include <math.h>
#include <conio.h>
double ed_f(double x,double y) { return(-y*y); }
double ed_solution(double x) { return(1.0/(1.0+x)); }
long ed_nystrom(double a,double b,double y,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 a,b,y0; double res[101][2]; clrscr(); nbpts=10; a=0;b=1;y0=1; printf("Méthode de Nystr?m\n"); for(n=10;n<=10000;n*=10) { nmax=ed_nystrom(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*n); } }

نتيجة البرنامج

  • في كل خطوة البرنامج يظهر:

  • قيمة xn

  • قيمة yn

  • الارتياب (الفرق بين yالمحسوبة والقيمة المضبوطةy(xn) )

ويظهر أيضا، لكل قيمة h، نسبة التقارب، العلاقة بين الارتياب وh2 في نقطة اعتباطية من المجال (نختار x = 1)

h =  0.1000
     xn            yn              الارتياب
  0.0000  1.0000000000000000e+00  0.00e+00
  0.1000  9.0975000000000006e-01  6.59e-04
  0.2000  8.3447098749999993e-01  1.14e-03
  0.3000  7.7048163420415505e-01  1.25e-03
  0.4000  7.1574259777081883e-01  1.46e-03
  0.5000  6.6802414095141105e-01  1.36e-03
  0.6000  6.2649134719204469e-01  1.49e-03
  0.7000  5.8952585933011048e-01  1.29e-03
  0.8000  5.5698319942826369e-01  1.43e-03
  0.9000  5.2747980244104153e-01  1.16e-03
  1.0000  5.0133621103161563e-01  1.34e-03
نسبة = 1.34e-01
h =  0.0100
     xn            yn              الارتياب
  0.0000  1.0000000000000000e+00  0.00e+00
  0.1000  9.0909846304181841e-01  7.55e-06
  0.2000  8.3334498858948347e-01  1.17e-05
  0.3000  7.6924454563473610e-01  1.38e-05
  0.4000  7.1430045347788507e-01  1.47e-05
  0.5000  6.6668168459206234e-01  1.50e-05
  0.6000  6.2501489371619001e-01  1.49e-05
  0.7000  5.8824983072805048e-01  1.45e-05
  0.8000  5.5556960654551146e-01  1.41e-05
  0.9000  5.2632929097502856e-01  1.35e-05
  1.0000  5.0001292828105715e-01  1.29e-05
     نسبة = 1.29e-01
h =  0.0010
     xn            yn              الارتياب
  0.0000  1.0000000000000000e+00  0.00e+00
  0.1000  9.0909098426961954e-01  7.52e-08
  0.2000  8.3333344916606789e-01  1.16e-07
  0.3000  7.6923090591630972e-01  1.37e-07
  0.4000  7.1428586023760565e-01  1.46e-07
  0.5000  6.6666681503825242e-01  1.48e-07
  0.6000  6.2500014675295934e-01  1.47e-07
  0.7000  5.8823543691185354e-01  1.43e-07
  0.8000  5.5555569309290365e-01  1.38e-07
  0.9000  5.2631592110111602e-01  1.32e-07
  1.0000  5.0000012546470518e-01  1.25e-07
     نسبة = 1.25e-01
h =  0.0001
     xn            yn              الارتياب
  0.0000  1.0000000000000000e+00  0.00e+00
  0.1000  9.0909090984227092e-01  7.51e-10
  0.2000  8.3333333449083291e-01  1.16e-09
  0.3000  7.6923077059640355e-01  1.37e-09
  0.4000  7.1428571574361899e-01  1.46e-09
  0.5000  6.6666666814837083e-01  1.48e-09
  0.6000  6.2500000146511059e-01  1.47e-09
  0.7000  5.8823529554275089e-01  1.43e-09
  0.8000  5.5555555692765746e-01  1.37e-09
  0.9000  5.2631579078623758e-01  1.31e-09
  1.0000  5.0000000125045985e-01  1.25e-09
  نسبة = 1.25e-01

الدالة الرئيسية ed_nystrom

معايير هذه الدالة هي كالتالي:

  • a و b : حدي المجال [ab]

  • y : الشرط البدئي y0 = y(a) = y(x0)i

  • n : عدد المجالات الصغيرة (ذات الطول h)

  • res : جدول لحفظ القيم (xnyn)

تستدعي الدالة الرئيسية الدالة التالية التي يجب تعريفها في البرنامج الرئيسي أيضا:

 double ed_f(double x,double y)

 تمثل هذه الدالة الدالة الرياضية  f .

يتم حساب القيم التقريبية لـ y(x)i من خلال طريقة نيستروم. وترجع بعدد القيم (xnyn) المحفوظة (أقصى حد هو 100).

سيتم حفظ هذه القيم في الجدول res : حيث سيحتوي العمود 0 على قيم xn بينما سيمثل العمود 1 قيم yn المقابلة.

ملاحظة: تتم الخطوة الأولى بمساعدة طريقة أولير المطورة.

long ed_nystrom(double a,double b,double y,long n,double res[101][2]) { double x,h,tol,xx,exx,yy,z; long i,j; h=(b-a)/n; tol=h/2; x=a; res[0][0]=a; res[0][1]=y; j=1; z=y; for(i=1;i<=n;i++) { if(i==1) y += h*ed_f(x+h/2,y+h/2*ed_f(x,y)); else { yy=z+2*h*ed_f(x,y); z=y; y=yy; } 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); } 

 دالة الإظهار ed_affichage

معايير دالة الإظهار هي كالتالي:

  • a و b : حدي المجال [ab]
    n : عدد المجالات الصغيرة (ذات الطول h)
    nmax : العدد الأقصى للقيم المحتفظ بها
    res : جدول يحفظ القيم (xnyn)
    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]));
}

}
}
]]>
aj.abderahim@gmail.com (أجديك عبد الرحيم) التكامل الرقمي للمعادلات التفاضلية العادية Tue, 18 Jun 2013 13:31:03 +0000
طريقة رونج كوطا ذات الرتبة 4 http://www.isnaha.com/isnaha_new/برمجها/item/950-طريقة-رونج-كوطا-ذات-الرتبة-4 http://www.isnaha.com/isnaha_new/برمجها/item/950-طريقة-رونج-كوطا-ذات-الرتبة-4
طريقة رونج كوطا ذات الرتبة 4
 Range-kutta
 

مصطلحات
العربية: طريقة رونج كوطا ذات الرتبة 4
الإنجليزية: Runge Kutte Method
الفرنسية: Méthode de Runge Kutta d'ordre 4

المبدأ

نجزئ المجال الاندماج [a, b] إلى N مجال ذي طول يساوي h = ،
محدودة بنقط التجزيء xn = a + nh , n = 0,... N .

لكل أفصول xn، نحسب القيمة المقربة yn للقيمة الصحيحة y(xn) لدالة y .

لهدا نبدل ميل المماس f(xn, yn) بمتوسط cefficientée لهدا الميل مع 3 قيمة مصححة على التوالي في 3 نقط من المجال [xnxn+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
k1 = f(x0, y0))
    = -y02
    = -12
    = -1
k2 = f(x0+(h/2), y0+h(k1/2))
    = -(y0+h(k1/2))2
    = -(1+(1/10)(-1/2))2
    = -361/400
k3 = f(x0+(h/2), y0+h(k2/2))
    = -(y0+h(k2/2))2
    = -(1+(1/10)(-361/800))2
    = -58354321/64000000
k4 = f(x0+h, y0+hk3)
    = -(y0+hk3)2
    = -(1+(1/10)(-58354321/64000000))2
    = -(581645679/640000000)2
y1 = y0 +(h/6)(k1+2k2+2k3+k4)
    = 1+(1/60)( -1-2(361/400)-2(58354321/64000000)-(581645679/640000000)2)


البرمجة

#include <math.h>
#include <conio.h>
double ed_f(double x,double y) { return(-y*y); }
double ed_solution(double x) { return(1.0/(1.0+x)); }
long ed_runge_kutta(double a,double b,double y,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 a,b,y0; double res[101][2]; clrscr(); nbpts=10; a=0;b=1;y0=1; printf("Méthode de Runge-Kutta\n"); for(n=10;n<=1000;n*=10) { nmax=ed_runge_kutta(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]))*pow(n,4)); } }

نتيجة البرنامج

  • في كل خطوة البرنامج يظهر:

  • قيمة xn

  • قيمة yn

  • الارتياب (الفرق بين yالمحسوبة والقيمة المضبوطةy(xn) )

ويظهر أيضا، لكل قيمة h، نسبة التقارب، العلاقة بين الارتياب و h4 في نقطة اعتباطية من المجال (نختار x = 1 )

h =  0.1000
     xn              yn            الارتياب
  0.0000  1.0000000000000000e+00  0.00e+00
  0.1000  9.0909118633221964e-01  2.77e-07
  0.2000  8.3333372884307211e-01  3.96e-07
  0.3000  7.6923120575328652e-01  4.37e-07
  0.4000  7.1428615389276151e-01  4.40e-07
  0.5000  6.6666709106586264e-01  4.24e-07
  0.6000  6.2500040094917386e-01  4.01e-07
  0.7000  5.8823566858083909e-01  3.74e-07
  0.8000  5.5555590318321424e-01  3.48e-07
  0.9000  5.2631611126425815e-01  3.22e-07
  1.0000  5.0000029758023101e-01  2.98e-07
   نسبة = 2.98e-03 
h =  0.0100
     xn              yn            الارتياب
  0.0000  1.0000000000000000e+00  0.00e+00
  0.1000  9.0909090911944701e-01  2.85e-11
  0.2000  8.3333333337395810e-01  4.06e-11
  0.3000  7.6923076927553513e-01  4.48e-11
  0.4000  7.1428571433074217e-01  4.50e-11
  0.5000  6.6666666671009622e-01  4.34e-11
  0.6000  6.2500000004099976e-01  4.10e-11
  0.7000  5.8823529415591624e-01  3.83e-11
  0.8000  5.5555555559106540e-01  3.55e-11
  0.9000  5.2631578950654223e-01  3.29e-11
  1.0000  5.0000000003037637e-01  3.04e-11
   نسبة= 3.04e-03
h =  0.0010
     xn              yn            الارتياب
  0.0000  1.0000000000000000e+00  0.00e+00
  0.1000  9.0909090909091184e-01  2.81e-15
  0.2000  8.3333333333333726e-01  4.03e-15
  0.3000  7.6923076923077394e-01  4.83e-15
  0.4000  7.1428571428571919e-01  5.05e-15
  0.5000  6.6666666666667118e-01  4.66e-15
  0.6000  6.2500000000000422e-01  4.38e-15
  0.7000  5.8823529411765119e-01  4.30e-15
  0.8000  5.5555555555555947e-01  4.10e-15
  0.9000  5.2631578947368829e-01  4.27e-15
  1.0000  5.0000000000000366e-01  3.83e-15
   نسبة = 3.83e-03

الدالة الرئيسية ed_runge_kutta

معايير هذه الدالة هي كالتالي:

  • a و b : حدي المجال [ab]

  • y : الشرط البدئي y0 = y(a) = y(x0)i

  • n : عدد المجالات الصغيرة (ذات الطول h)

  • res : جدول لحفظ القيم (xnyn)

تستدعي الدالة الرئيسية الدالة التالية التي يجب تعريفها في البرنامج الرئيسي أيضا:

 double ed_f(double x,double y)

 تمثل هذه الدالة الدالة الرياضية  f .

يتم حساب القيم التقريبية لـ y(x)i من خلال طريقة رونج كيتا. وترجع بعدد القيم (xnyn) المحفوظة (أقصى حد هو 100).

سيتم حفظ هذه القيم في الجدول res :حيث سيحتوي العمود 0 على قيم xn بينما سيمثل العمود 1 قيم yn المقابلة.

long ed_runge_kutta(double a,double b,double y,long n,double res[101][2]) 
{
double x,h,k1,k2,k3,k4,tol,xx,exx;
long i,j;
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++)
{
k1=ed_f(x,y);
k2=ed_f(x+h/2,y+h*k1/2);
k3=ed_f(x+h/2,y+h*k2/2);
k4=ed_f(x+h,y+h*k3);
y += h*(k1+2*k2+2*k3+k4)/6;
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);
}

دالة الإظهار ed_affichage


معايير دالة الإظهار هي كالتالي:

  • a و b : حدي المجال [ab]
    n : عدد المجالات الصغيرة (ذات الطول h)
    nmax : العدد الأقصى للقيم المحتفظ بها
    res : جدول يحفظ القيم (xnyn)
    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]));
}

}
}
]]>
aj.abderahim@gmail.com (أجديك عبد الرحيم) التكامل الرقمي للمعادلات التفاضلية العادية Tue, 18 Jun 2013 13:30:11 +0000
طريقة أولير كوشي http://www.isnaha.com/isnaha_new/برمجها/item/949-طريقة-أولير-كوشي http://www.isnaha.com/isnaha_new/برمجها/item/949-طريقة-أولير-كوشي
طريقة أولير كوشي
 Euler-Cauchy
 

مصطلحات
العربية: طريقة أولير كوشي
الإنجليزية: Euler Cauchy Method
الفرنسية: Méthode d'Euler Cauchy

المبدأ

حل رقميا  y'(x) = f(xy(x))  

بحيث

 x  [ab]   y(a) = y0

نجزئ مجال الاندماج [a, b] إلى N مجال ذي طول يساوي h = ،
محدودة بنقط التجزيء xn = a + nh , n = 0,... N .

لكل أفصول xn ،نحسب القيمة المقربةyn  للقيمة الصحيحةy(xn) لدالة y.

لهدا نبدل ميل المماس f(xn, yn) بمتوسط هدا الميل مع القيمة المصححة فيxn+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/2)( f(x0, y0)+f(x0+h, y0+h f(x0, y0)))
     = y0 - (h/2)(y02+( y0 -hy02 )2)
     = 1 - (1/20)(12+(1-(1/10)12)2)
     = 1819/2000

x2 = 2/10
y2 = y1 + (h/2)( f(x1, y1)+f(x1+h, y1+h f(x1, y1)))
     = y1 - (h/2)(y12+( y1 -hy12 )2)
     = 1819/2000 - (1/20)((1819/2000)2+(1819/2000-(1/10)(1819/2000)2)2)


البرمجة

#include <math.h>
#include <conio.h>

double ed_f(double x,double y) { return(-y*y); }
double ed_solution(double x) { return(1.0/(1.0+x)); }
long ed_euler_cauchy(double a,double b,double y,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 a,b,y0; double res[101][2]; clrscr(); nbpts=10; a=0;b=1;y0=1; printf("Méthode d'Euler-Cauchy\n"); for(n=10;n<=10000;n*=10) { nmax=ed_euler_cauchy(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*n); } }

نتيجة البرنامج

  • في كل خطوة البرنامج يظهر:

  • قيمة xn

  • قيمة yn

  • الارتياب (الفرق بين ynالمحسوبة والقيمة المضبوطةy(xn) )

ويظهر أيضا ،لكل قيمة h، نسبة التقارب، العلاقة بين الارتياب وh2 في نقطة اعتباطية من المجال (نختار x = 1 )

h =  0.1000
      xn              yn           الارتياب  
  0.0000  1.0000000000000000e+00  0.00e+00
  0.1000  9.0949999999999998e-01  4.09e-04
  0.2000  8.3396214846890249e-01  6.29e-04
  0.3000  7.6997115403838756e-01  7.40e-04
  0.4000  7.1507467430115557e-01  7.89e-04
  0.5000  6.6746716934952444e-01  8.01e-04
  0.6000  6.2579033562771691e-01  7.90e-04
  0.7000  5.8900298033961918e-01  7.68e-04
  0.8000  5.5629374665503883e-01  7.38e-04
  0.9000  5.2702111163636556e-01  7.05e-04
  1.0000  5.0067122128275432e-01  6.71e-04
      نسبة = 6.71e-02
h =  0.0100
      xn              yn           الارتياب  
  0.0000  1.0000000000000000e+00  0.00e+00
  0.1000  9.0909470137414994e-01  3.79e-06
  0.2000  8.3333917324553641e-01  5.84e-06
  0.3000  7.6923765696640956e-01  6.89e-06
  0.4000  7.1429306526662695e-01  7.35e-06
  0.5000  6.6667413571001111e-01  7.47e-06
  0.6000  6.2500738366174291e-01  7.38e-06
  0.7000  5.8824247460205625e-01  7.18e-06
  0.8000  5.5556246758329719e-01  6.91e-06
  0.9000  5.2632240025000077e-01  6.61e-06
  1.0000  5.0000629686925135e-01  6.30e-06
      نسبة = 6.30e-02
h =  0.0010
      xn              yn           الارتياب  
  0.0000  1.0000000000000000e+00  0.00e+00
  0.1000  9.0909094669249346e-01  3.76e-08
  0.2000  8.3333339125673489e-01  5.79e-08
  0.3000  7.6923083756607158e-01  6.83e-08
  0.4000  7.1428578723447389e-01  7.29e-08
  0.5000  6.6666674080246002e-01  7.41e-08
  0.6000  6.2500007330169050e-01  7.33e-08
  0.7000  5.8823536541378374e-01  7.13e-08
  0.8000  5.5555562419600368e-01  6.86e-08
  0.9000  5.2631585513098234e-01  6.57e-08
  1.0000  5.0000006254687424e-01  6.25e-08
      نسبة = 6.25e-02
h =  0.0001
      xn              yn           الارتياب  
  0.0000  1.0000000000000000e+00  0.00e+00
  0.1000  9.0909090946660254e-01  3.76e-10
  0.2000  8.3333333391208975e-01  5.79e-10
  0.3000  7.6923076991357764e-01  6.83e-10
  0.4000  7.1428571501463878e-01  7.29e-10
  0.5000  6.6666666740746761e-01  7.41e-10
  0.6000  6.2500000073247952e-01  7.32e-10
  0.7000  5.8823529483009662e-01  7.12e-10
  0.8000  5.5555555624147845e-01  6.86e-10
  0.9000  5.2631579012980312e-01  6.56e-10
  1.0000  5.0000000062504335e-01  6.25e-10
      نسبة = 6.25e-02

الدالة الرئيسية ed_euler_cauchy
 
معايير هذه الدالة هي كالتالي:
 a و b : حدي المجال [ab]
y : الشرط البدئي y0 = y(a) = y(x0)i
n : عدد المجالات الصغيرة (ذات الطول h)
res : جدول لحفظ القيم (xnyn)
 
 تستدعي الدالة الرئيسية الدالة التالية التي يجب تعريفها في البرنامج الرئيسي أيضا:
  double ed_f(double x,double y)
تمثل هذه الدالة الدالة الرياضية  f.
 يتم حساب القيم التقريبية لـ y(x)i من خلال طريقة أولير كوشي. وترجع بعدد القيم (xnyn) المحفوظة (أقصى حد هو 100). سيتم حفظ هذه القيم في الجدول res :حيث سيحتوي العمود 0 على قيم xn بينما سيمثل العمود 1 قيم yn المقابلة. 
long ed_euler_cauchy(double a,double b,double y,long n,double res[101][2])
{
double x,h,fxy,tol,xx,exx;
long i,j;
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++)
{
fxy=ed_f(x,y);
y += h*(fxy+ed_f(x+h,y+h*fxy))/2;
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);
}

 دالة الإظهار ed_affichage

معايير دالة الإظهار هي كالتالي:

  • a و b : حدي المجال [ab]
    n : عدد المجالات الصغيرة (ذات الطول h)
    nmax : العدد الأقصى للقيم المحتفظ بها
    res : جدول يحفظ القيم (xnyn)
    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]));
}

}
}
]]>
aj.abderahim@gmail.com (أجديك عبد الرحيم) التكامل الرقمي للمعادلات التفاضلية العادية Tue, 18 Jun 2013 13:29:14 +0000
طريقة أولير المطورة http://www.isnaha.com/isnaha_new/برمجها/item/948-طريقة-أولير-المطورة http://www.isnaha.com/isnaha_new/برمجها/item/948-طريقة-أولير-المطورة
طريقة أولير المطورة
Euler method 
 

مصطلحات
العربية: طريقة أولير المطورة
الإنجليزية: Improved Euler's Method
الفرنسية: Méthode d'Euler ameliorée

المبدأ

نجزئ المجال الاندماج [a, b] إلى N مجال ذي طول يساوي h = ،
محدد بنقط التجزيء xn = a + nh , n = 0,... N .

لكل أفصول xn ،نحسب القيمة المقربة yn  للقيمة الصحيحة y(xn) لدالة y .

لهدا نبدل ميل المماس f(xn, yn) بـ القيمة المصححة في وسط المجال [xn, xn+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/2), y0+(h/2) f(x0, y0))
    = y0 - h( y0 -(h/2)y02 )2
    = 1 - (1/10)(1-(1/20)12)2
    = 3639/4000

x2 = 2/10
y2 = y1 + h f(x1+(h/2), y1+(h/2) f(x1, y1))
    = y1 - h ( y1 -(h/2)y12 )2
    = 3639/4000 - (1/10)(3639/4000-(1/20)(3639/4000)2)2


البرمجة

#include <math.h>
#include <conio.h>

double ed_f(double x,double y)
{
 return(-y*y);
}

double ed_solution(double x)
{
 return(1.0/(1.0+x));
}

long ed_tangente_amelioree(double a,double b,double y,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 a,b,y0;
 double res[101][2];
 clrscr();
 nbpts=10;
 a=0;b=1;y0=1;
 printf("Méthode d'Euler améliorée\n");
 for(n=10;n<=10000;n*=10)
 {
  nmax=ed_tangente_amelioree(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*n);
 }
}

نتيجة البرنامج

  • في كل خطوة البرنامج يظهر:

  • قيمة xn

  • قيمة yn

  • الارتياب (الفرق بين yالمحسوبة والقيمة المضبوطةy(xn) )

ويظهر أيضا ،لكل قيمة h، نسبة التقارب، العلاقة بين الارتياب و h2 في نقطة اعتباطية من المجال (نختارx = 1 )

h =  0.1000
    xn              yn            الارتياب
  0.0000  1.0000000000000000e+00  0.00e+00
  0.1000  9.0975000000000006e-01  6.59e-04
  0.2000  8.3434374561872360e-01  1.01e-03
  0.3000  7.7041776079476987e-01  1.19e-03
  0.4000  7.1554809944410103e-01  1.26e-03
  0.5000  6.6794532411677843e-01  1.28e-03
  0.6000  6.2626051029611207e-01  1.26e-03
  0.7000  5.8945803964048626e-01  1.22e-03
  0.8000  5.5672991484754586e-01  1.17e-03
  0.9000  5.2743665240683268e-01  1.12e-03
  1.0000  5.0106563581429009e-01  1.07e-03
     نسبة = 1.07e-01
h =  0.0100
     xn              yn           الارتياب
  0.0000  1.0000000000000000e+00  0.00e+00
  0.1000  9.0909663426283682e-01  5.73e-06
  0.2000  8.3334214754473601e-01  8.81e-06
  0.3000  7.6924116268653298e-01  1.04e-05
  0.4000  7.1429680471874069e-01  1.11e-05
  0.5000  6.6667793341463788e-01  1.13e-05
  0.6000  6.2501113639096484e-01  1.11e-05
  0.7000  5.8824612272614063e-01  1.08e-05
  0.8000  5.5556597816773001e-01  1.04e-05
  0.9000  5.2632575684848548e-01  9.97e-06
  1.0000  5.0000949324076738e-01  9.49e-06
     نسبة = 9.49e-02
h =  0.0010
     xn              yn            الارتياب
  0.0000  1.0000000000000000e+00  0.00e+00
  0.1000  9.0909096552923185e-01  5.64e-08
  0.2000  8.3333342027161306e-01  8.69e-08
  0.3000  7.6923087179426586e-01  1.03e-07
  0.4000  7.1428582377147642e-01  1.09e-07
  0.5000  6.6666677793223084e-01  1.11e-07
  0.6000  6.2500011001218381e-01  1.10e-07
  0.7000  5.8823540111855599e-01  1.07e-07
  0.8000  5.5555565856969691e-01  1.03e-07
  0.9000  5.2631588800981521e-01  9.85e-08
  1.0000  5.0000009386729316e-01  9.39e-08
    نسبة = 9.39e-02
h =  0.0001
     xn              yn            الارتياب
  0.0000  1.0000000000000000e+00  0.00e+00
  0.1000  9.0909090965448325e-01  5.64e-10
  0.2000  8.3333333420152000e-01  8.68e-10
  0.3000  7.6923077025504372e-01  1.02e-09
  0.4000  7.1428571537916319e-01  1.09e-09
  0.5000  6.6666666777792982e-01  1.11e-09
  0.6000  6.2500000109877951e-01  1.10e-09
  0.7000  5.8823529518638029e-01  1.07e-09
  0.8000  5.5555555658449340e-01  1.03e-09
  0.9000  5.2631579045791688e-01  9.84e-10
  1.0000  5.0000000093761643e-01  9.38e-10
    نسبة = 9.38e-02

الدالة الرئيسية ed_tangente_amelioree

معايير هذه الدالة هي كالتالي:

 a و b : حدي المجال [ab]
y : الشرط البدئي y0 = y(a) = y(x0)i
n : عدد المجالات الصغيرة (ذات الطول h)
res : جدول لحفظ القيم (xnyn)

تستدعي الدالة الرئيسية الدالة التالية التي يجب تعريفها في البرنامج الرئيسي أيضا:

  double ed_f(double x,double y)

تمثل هذه الدالة الدالة الرياضية  f.

 يتم حساب القيم التقريبية لـ y(x)i من خلال طريقة أولير. وترجع بعدد القيم (xnyn) المحفوظة (أقصى حد هو 100).

 سيتم حفظ هذه القيم في الجدول res :حيث سيحتوي العمود 0 على قيم xn بينما سيمثل العمود 1 قيم yn المقابلة. 

long ed_tangente_amelioree(double a,double b,double y,long n,double res[101][2]) 
{ double x,h,tol,xx,exx;
long i,j;
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++)
{
y += h*ed_f(x+h/2,y+h/2*ed_f(x,y));
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);
}

  • دالة الإظهار 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]));
}

}
}
]]>
aj.abderahim@gmail.com (أجديك عبد الرحيم) التكامل الرقمي للمعادلات التفاضلية العادية Tue, 18 Jun 2013 13:27:32 +0000
طريقة أولير http://www.isnaha.com/isnaha_new/برمجها/item/947-طريقة-أولير http://www.isnaha.com/isnaha_new/برمجها/item/947-طريقة-أولير
طريقة أولير
 Euler method
 

مصطلحات
العربية: طريقة أولير
الإنجليزية: Euler Method
الفرنسية: Méthode d'Euler

المبدأ

نجزئ المجال [a, b] إلى N مجال مصغر ذي طول يساوي h = ،
المجالات المصغرة محدودة بنقط التجزيء xn = a + nh , n = 0,... N.

لكل أفصول xn، نحسب القيمة المقربة yn للقيمة الصحيحة y(xn) لدالة y.

لهذا، نعتبر مماس في (xn, yn) للمنحنى التام المار من نقطه، القيمةyn+1 هي أرتوب نقطة هدا المستقيم ذي أفصول xn+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, y0)
     = y0 - h y02
     = 1 - (1/10)12
     = 9/10

x2 = 2/10
y2 = y1 + h f(x1, y1)
     = y1 - h y12
     = (9/10) - (1/10)(9/10)2
     = 819/1000


البرمجة

#include <math.h>
#include <conio.h>

double ed_f(double x,double y) { return(-y*y); }
double ed_solution(double x) { return(1.0/(1.0+x)); }
long ed_euler(double a,double b,double y,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 a,b,y0; double res[101][2]; clrscr(); nbpts=10; a=0;b=1;y0=1; printf("Méthode d'Euler\n"); for(n=10;n<=10000;n*=10) { nmax=ed_euler(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
    الارتياب (الفرق بين yالمحسوبة والقيمة المضبوطةy(xn) )

ويظهر أيضا، لكل قيمة، نسبة التقارب، العلاقة بين الارتياب وفي نقطة إعتباطية من المجال (نختارx = 1 )

    xn              yn            الارتياب
  0.0000  1.0000000000000000e+00  0.00e+00
  0.1000  9.0000000000000002e-01  9.09e-03
  0.2000  8.1900000000000006e-01  1.43e-02
  0.3000  7.5192390000000009e-01  1.73e-02
  0.4000  6.9538494486087910e-01  1.89e-02
  0.5000  6.4702892270696233e-01  1.96e-02
  0.6000  6.0516428002502909e-01  1.98e-02
  0.7000  5.6854189944320788e-01  1.97e-02
  0.8000  5.3621791030095878e-01  1.93e-02
  0.9000  5.0746494556820609e-01  1.89e-02
  1.0000  4.8171287847015187e-01  1.83e-02
      نسبة = 1.83e-01
h =  0.0100
     xn              yn            الارتياب
  0.0000  1.0000000000000000e+00  0.00e+00
  0.1000  9.0829282588834559e-01  7.98e-04
  0.2000  8.3205259796437259e-01  1.28e-03
  0.3000  7.6766248975763363e-01  1.57e-03
  0.4000  7.1255342431233093e-01  1.73e-03
  0.5000  6.6484990411095146e-01  1.82e-03
  0.6000  6.2315052002916993e-01  1.85e-03
  0.7000  5.8638691400614118e-01  1.85e-03
  0.8000  5.5373029959869036e-01  1.83e-03
  0.9000  5.2452781589306963e-01  1.79e-03
  1.0000  4.9825816164586700e-01  1.74e-03
      نسبة = 1.74e-01
h =  0.0010
     xn              yn            الارتياب
  0.0000  1.0000000000000000e+00  0.00e+00
  0.1000  9.0901203784010098e-01  7.89e-05
  0.2000  8.3320657680096466e-01  1.27e-04
  0.3000  7.6907536755137540e-01  1.55e-04
  0.4000  7.1411389036029838e-01  1.72e-04
  0.5000  6.6648631430783711e-01  1.80e-04
  0.6000  6.2481627065328638e-01  1.84e-04
  0.7000  5.8805156369036982e-01  1.84e-04
  0.8000  5.5537402970003502e-01  1.82e-04
  0.9000  5.2613789141912859e-01  1.78e-04
  1.0000  4.9982662405687811e-01  1.73e-04
      نسبة = 1.73e-01
h =  0.0001
     xn              yn            الارتياب
  0.0000  1.0000000000000000e+00  0.00e+00
  0.1000  9.0908303119289890e-01  7.88e-06
  0.2000  8.3332067067247628e-01  1.27e-05
  0.3000  7.6921524315425815e-01  1.55e-05
  0.4000  7.1426854579220944e-01  1.72e-05
  0.5000  6.6664864454004780e-01  1.80e-05
  0.6000  6.2498163914253235e-01  1.84e-05
  0.7000  5.8821693205906611e-01  1.84e-05
  0.8000  5.5553741289040837e-01  1.81e-05
  0.9000  5.2629800859634635e-01  1.78e-05
  1.0000  4.9998267042955358e-01  1.73e-05
     نسبة = 1.73e-01
 

الدالة الرئيسية ed_euler

معايير هذه الدالة هي كالتالي:

  • a و b : حدي المجال [ab]

  • y : الشرط البدئي y0 = y(a) = y(x0)i

  • n : عدد المجالات الصغيرة (ذات الطول h)

  • res : جدول لحفظ القيم (xnyn)

تستدعي الدالة الرئيسية الدالة التالية التي يجب تعريفها في البرنامج الرئيسي أيضا:

 double ed_f(double x,double y)

 تمثل هذه الدالة الدالة الرياضية  f .

يتم حساب القيم التقريبية لـ y(x)i من خلال طريقة أولير. وترجع بعدد القيم (xnyn) المحفوظة (أقصى حد هو 100) ;

سيتم حفظ هذه القيم في الجدول res : حيث سيحتوي العمود 0 على قيم xn بينما سيمثل العمود 1 قيم yn المقابلة.

long ed_euler(double a,double b,double y,long n,double res[101][2])
{
double x,h,tol,xx,exx;
long i,j;

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++)
{
y += h*ed_f(x,y);
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);
}

 دالة الإظهار ed_affichage

معايير دالة الإظهار هي كالتالي:

  • a و b : حدي المجال [ab]
    n : عدد المجالات الصغيرة (ذات الطول h)
    nmax : العدد الأقصى للقيم المحتفظ بها
    res : جدول يحفظ القيم (xnyn)
    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]));
}

}
}
]]>
aj.abderahim@gmail.com (أجديك عبد الرحيم) التكامل الرقمي للمعادلات التفاضلية العادية Tue, 18 Jun 2013 13:26:33 +0000
التكامل الرقمي للمعادلات التفاضلية العادية http://www.isnaha.com/isnaha_new/برمجها/item/946-التكامل-الرقمي-للمعادلات-التفاضلية-العادية http://www.isnaha.com/isnaha_new/برمجها/item/946-التكامل-الرقمي-للمعادلات-التفاضلية-العادية
التكامل الرقمي للمعادلات التفاضلية العادية 
 
 

الشروط البدئية
 
·  ليكن [a,bمجالا مغلقا، التطبيق ،f :  هو تطبيق متفاضل .y : 
·  نعتبر المسألة التالية (مسألة كوشي): x[ab] وy'(x) = f(xy(x))  مع y(a)=y0 (الحالة البدئية).
 (f دالة معطاة وy دالة مجهولة مراد تحديدها ،  y معطاة)

 الطرق ذات الخطوات المتفرقة 

الطرق ذات الخطوات المرتبطة 

طريقة أوليرEuler 

طريقة أولير المتطورة Euler amélioré 

طريقة أولير  -كوشي Euler-Cauchy 

طريقة رونج-كيتا Runge-Kutta ذات الرتبة 4

طريقة نيستروم Nystr?m   

طريقة أولير الضمنية Euler implicite 

طريقة آدم-باشفورة-مولت Adams-Bashforth-Moulton

. 


مشكلة الحالة في المحدات

 

ليكن  [a, b]مجال مغلق، التطبيق   بحيث f :   هو تطبيق متفاضل مرتين   .y :

نعتبر المسألة التالية: x[a, b] و y''(x) = f(x) بحيثy(a)=y0   و y(b)=yN  (الحالة في المحدات)

(f دالة معطاة،y   دالة مجهولة مراد تحديدها  تنتمتي إلى y0 و yN معطاة)

   الطريقة الرقمية 

طريقة الفوارق المنتهية 

]]>
aj.abderahim@gmail.com (أجديك عبد الرحيم) التكامل الرقمي للمعادلات التفاضلية العادية Tue, 18 Jun 2013 13:25:57 +0000