โปรแกรมช่วยหาคำตอบระบบสมการไม่เชิงเส้นด้วย Newton-Raphson Method

สวัสดีครับ ในบทความนี้ จะขอนำเสนอการประยุกต์คอนสตรัคเตอร์แบบที่ใช้การโหลดสมาชิกของเมตริกซ์ จากไดนามิกส์อาเรย์  โดยจะนำเสนอผ่านระเบียบวิธีการหาคำตอบของระบบสมการไม่เชิงเส้นด้วย Newton-Raphson Method ผมจะขออธิบายแบบสั้นๆก่อนนำเข้าสู่โค้ดโปรแกรมภาษา C++ โดยปกติระบบสมการเชิงเส้น เราสามารถเขียนอยู่ในรูปเมตริกซ์ได้เป็น

[A][x] = [b]

โดย เมตริกซ์ A เป็นเมตริกซ์ที่มีสมาชิกของเมตริกซ์เป็นค่าคงที่ ดังแสดงในรูปที่ 1

ซึ่งเมตริกซ์ A จะเป็นเมตริกซ์ที่คงที่ ดังนั้นจึงสามารถใช้คลาส CMatrix อ่านข้อมูลสมาชิกของเมตริกซ์ A และ b ดังแสดงในบทความ การหาคำตอบของระบบสมการเชิงเส้นได้ แต่ในกรณีเป็นระบบสมการไม่เชิงเส้น (Nonlinear Equation System) สมาชิกของเมตริกซ์ A จะเป็นฟังก์ชันของตัวแปรที่ต้องการหาคำตอบ ดังแสดงในรูปที่ 2 การหาคำตอบของระบบสมการไม่เชิงเส้นจึงต้องใช้การหาคำตอบของระบบสมการด้วย Newton-Raphson Method 

Newton-Raphson Method จะมีขั้นตอนดำเนินการดังนี้

1. เขียนสมการทั้งหมดให้อยู่ในรูปแบบดังนี้

จากตัวอย่าง 3 สมการจะสามารถเขียน f1 , f2 , f3 ได้ดังนี้

2. การหาคำตอบของระบบสมการ จะเริ่มจากการเขียนระบบสมการในรูปแบบของอนุกรมเทย์เลอร์ โดยสามารถเขียนในรูปเมตริกซ์ได้ดังนี้

               

เมตริกซ์ J จะเรียกว่า Jacobian matrix ซึ่งสมาชิกต่างๆหาได้ดังนี้

จากตัวอย่าง 3 สมการ จะเขียนระบบสมการได้ดังนี้

3. ทำการหาคำตอบของระบบสมการนี้ โดยการกำหนดค่าเริ่มต้น  

แทนลงในเมตริกซ์ J และ f และหาคำตอบของระบบสมการในข้อ 2 ดังกล่าว

4. คำนวณหาผลลัพธ์ใหม่จาก  

5. ตรวจสอบผลลัพธ์โดยการแทนคำตอบจากข้อ 4 ลงใน เมตริกซ์ f และตรวจสอบค่าของสมาชิกในเมตริกซ์ โดยหากคำตอบที่แทนลงในเมตริกซ์ f ทำให้สมาชิกในเมตริกซ์ f มีค่าเท่ากับ 0  จะได้ว่า คำตอบดังกล่าวเป็นคำตอบของระบบสมการ แต่ในทางปฏิบัติจะกำหนดเงื่อนไขการลู่เข้าหาคำตอบไว้ โดยหากค่าในเมตริกซ์ f ทุกค่าต่ำกว่าที่กำหนดจะถือว่าชุดคำตอบนี้เป็นคำตอบของระบบสมการ แต่หากไม่ใช่จะย้อนกลับไปหาคำตอบของระบบสมการในข้อ 2 ต่อไปและจะทำซ้ำจนกระทั่งลู่เข้าสู่คำตอบของระบบสมการ

จากที่ได้อธิบายมาคร่าวๆ เราจะเห็นจุดที่จะนำเอาสมาชิกเมตริกซ์จากอาเรย์ 2D เข้ามาประยุกใช้ได้ นั่นคือขั้นตอนที่ 2 จากนั้นกำหนดคำตอบเริ่มต้นในข้อ 3 และทำการหาคำตอบของระบบสมการตามข้อ 2 โดยใช้ฟังก์ชัน SolveLU หรือ SolveGauss ได้ เรามาดู โค้ดโปรแกรมกันเลยครับ

#include "stdafx.h"

#include "Matrix.h"

#include "math.h"

int _tmain(int argc, _TCHAR* argv[])

{

        int n ;

int count ;

int row , col ;

double **J ;

double **f ;

double **x ;

double tol = 0.000001 ;

row = 3 ;

col = 3 ;

count = 0 ;

J = new double*[row];

f = new double*[row];

x = new double*[row];

for(int i = 0 ; i < row ; i++)

{

J[i] = new double[col] ;

f[i] = new double[col] ;

x[i] = new double[col] ;

}

//เดาคำตอบเริ่มต้น

x[0][0] = 1 ;   //x

x[1][0] = 1 ; //y

x[2][0] = 1 ; //z

       //คำนวณสมาชิกเมตริกซ f เริ่มต้น

f[0][0] = -(4*x[0][0]+x[1][0]*x[1][0] +x[2][0] -11) ;  //-f1:-(4x+y^2+z-11)

f[1][0] = -(x[0][0]+4*x[1][0] +x[2][0]*x[2][0] -18) ;  //-f2:-(x+4y+z^2-18)

f[2][0] = -(x[0][0]*x[0][0]+x[1][0] +4*x[2][0] -15) ;  //-f3:-(x^2+y+4z-15)

do

{      

               count++;

       //คำนวณเมตริซ์ J

      J[0][0] = 4 ;

      J[0][1] = 2*x[1][0] ;  //2y

      J[0][2] = 1 ;

      J[1][0] = 1 ;

      J[1][1] = 4 ;

      J[1][2] = 2*x[2][0] ; //2z

      J[2][0] = 2*x[0][0] ; //2x

      J[2][1] = 1 ;

      J[2][2] = 4 ;

     CMatrix A(3,3,J) ;

     CMatrix b(3,1,f) ;

     CMatrix dx ;                      //สร้าง object dx ไว้สำหรับเก็บคำตอบที่ Solve ได้

     dx.SolveLU(&A,&b);       // Solve ระบบสมการและเก็บคำตอบไว้ที่เมตริกซ dx

       

             // กำหนดคำตอบที่ได้ลงในเมตริกซ์ x 

     x[0][0] += dx.GetData(0,0);

     x[1][0] += dx.GetData(1,0);

     x[2][0] += dx.GetData(2,0);

             //คำนวณสมาชิกเมตริกซ f ใหม่

     f[0][0] = -(4*x[0][0]+x[1][0]*x[1][0] +x[2][0] -11) ;  //-f1:-(4x+y^2+z-11)

     f[1][0] = -(x[0][0]+4*x[1][0] +x[2][0]*x[2][0] -18) ;  //-f2:-(x+4y+z^2-18)

     f[2][0] = -(x[0][0]*x[0][0]+x[1][0] +4*x[2][0] -15) ;  //-f3:-(x^2+y+4z-15)     

}

        // check เงื่อนไขการลู่เข้าของคำตอบ

while (fabs(f[0][0]) > tol || fabs(f[1][0]) > tol || fabs(f[2][0]) > tol);   

printf("inter %d\tsolution are :\n",count);

printf("x = %.3f\n",x[0][0]);

printf("y = %.3f\n",x[1][0]);

printf("z = %.3f\n\n",x[2][0]);

printf("f1 = %.7f\n",f[0][0]);

printf("f2 = %.7f\n",f[1][0]);

printf("f3 = %.7f\n\n",f[2][0]);

// คืน memory ให้ระบบ

for (int i = 0 ; i < row ; i++)

{

delete [] J[i] ;

delete [] f[i] ;

delete [] x[i] ;

}

delete [] J ;

delete [] f ;

delete [] x ;

J = NULL ;

f = NULL ;

x = NULL ;

        

return 0 ;

}

ผลการประมวลผลโปรแกรมภาษา c++ ดังกล่าวสำหรับใช้หาคำตอบระบบสมการที่ไม่เป็นเชิงเส้นตามตัวอย่างแสดงได้ดังรูป

คำตอบของระบบสมการ

จากตัวอย่างโค้ดจะเห็นว่าในส่วนของการสร้างเมตริกซ์ J และ f เป็นการประยุกต์ใช้วิธีการประกาศคอนสตรัคเตอร์ตามที่ได้นำเสนอในบทความที่ผ่านมาครับ ซึ่งก็จะช่วยให้วิศวกรสามารถแก้ปัญหาในระบบต่างๆได้โดยง่าย

         อยากจะกล่าวให้เห็นปัญหาของการหาคำตอบของระบบสมการที่มากกว่า 3 สมการ ซึ่งท่านต้องสร้างเมตริกซ์ J และ f มากขึ้น ซึ่งการสร้างเมตริกซ์ตามที่ผมได้นำเสนออาจจะทำได้ช้าและยุ่งยาก ซึ่งเราอาจจะประยุกต์การสร้างฟังก์ชันสำหรับคำนวณ สมาชิกในเมตริกซ์ขึ้นมาซึ่งทำได้รวดเร็วและเป็นแบบแผน ทำให้เราสามารถประยุกต์ใช้โค้ดโปรแกรมเดียวกันนี้กับระบบสมการหลายๆระบบได้โดยการเปลี่ยนแปลงโค้ดเพียงเล็กน้อย ซึ่งจะข้อนำเสนอในบทความถัดๆไปนะครับ เพื่อให้เห็นประโยชน์ของการเขียนโปรแกรมแบบ OOP 

ตัวอย่างการโหลดสมาชิกเมตริกซ์จากอาเรย์ 2D มาใช้ใน Class CMatrix ในโปรแกรมภาษา C++

สวัสดีครับ บทความนี้ผมจะขอนำเสนอการโหลดสมาชิกของเมตริกซ์จากอาเรย์ 2D มาใช้ในคลาส CMatrix กันครับ เริ่มจาก จากบทความแนะนำ คลาส CMatrix  ได้อธิบาย คอนสตรัคเตอร์ตัวหนึ่งที่ใช้ในการโหลดสมาชิกของเมตริกซ์จากอาเรย์ 2D ซึ่งมีรูปแบบการประกาศ Object ดังนี้ครับ

CMatrix A(const int r , const int c , double **pD)

มาพิจารณาพารามิเตอร์ที่ต้องส่งผ่านคอนสตรัคเตอร์ตัวนี้ ซึ่งจะประกอบไปด้วย

1. จำนวนแถวของเมตริซ์ r

2. จำนวนหลักของเมตริกซ์ c

3. อาเรย์ 2D ซึ่งต้องเป็นแบบ pointer 

สมมุติตัวอย่าง ต้องการโหลดสมาชิกเมตริกซ์ขนาด 3x3 จากตัวแปรอาเรย์ชื่อ A เราจะเขียน code โปรแกรมภาษา C++ ในส่วนนี้ได้ดังนี้

        int row , col ;

double** data ;

row = 3 ;

col = 3 ;

data = new double*[row];

for (int i = 0 ; i < row ; i++)

{

data[i] = new double[col] ;

}

for (int j = 0 ; j < row ; j++)

{

for(int k = 0 ; k < col ; k++)

data[j][k] = j+k ;

}

cout<<"Matrix A"<<endl;

CMatrix A(row,col,data);

A.Show();

cout<<endl ;

cout<<"Matrix B"<<endl;

CMatrix B(3,3,data) ;

B.Show();

ผลลัพธ์การทำงานใน code ส่วนนี้จะแสดงได้ดังนี้

ผลการทำงานของโค้ดโปรแกรม

จากตัวอย่างโค้ดโปรแกรมภาษา C++ ที่ประมวลผลเมตริกซ์นี้ ผมได้นำเสนอ 2 รูปแบบ คือ ท่านสามารถผ่านขนาดของเมตริกซ์ด้วยตัวแปรก็ได้ หรือจะผ่านด้วยค่าก็ได้ เนื่องจาก คอนสตรัคเตอร์นี้ กำหนดชนิดของตัวแปรขนาดเป็นแบบ คงที่ (const int r , const int c , double **pD) เอาหล่ะครับ จบตรงนี้เราก็ได้ object ของคลาส CMatrix ที่พร้อมจะนำไปใช้ประมวลผลหรือดำเนินการทางเมตริกซ์ต่อ ไม่ว่าจะนำไป บวก ลบ หรือคูณเมตริกซ์ หรือแม้แต่จะนำไปหาคำตอบของระบบสมการดังที่เคยนำเสนอไปแล้ว ในบทความต่อไปจะนำ ตัวอย่างการใช้คอนสตรัคเตอร์ดังกล่าวในการแก้ปัญหาระบบสมการไม่เชิงเส้น โดยใช้ระเบียบวิธี Newton-Raphson 

โปรแกรมคอมพิวเตอร์ช่วยหาคำตอบของระบบสมการเชิงเส้น ด้วย Class CMatrix

สวัสดีครับ บทความ โปรแกรมภาษา C++ บทนี้จะขอนำเสนอการประยุกต์ใช้คลาส CMatrix ช่วยในการหาคำตอบของระบบสมการเชิงเส้น เราจะประดิษฐ์โปรแกรมสำหรับช่วยหาคำตอบของระบบสมการเชิงเส้น เพื่อให้ท่านผู้อ่่านนำไปใช้งาน ทั้ง นร นศ ที่ต้องการนำไปตรวจเช็คคำตอบกับวิธีการหาคำตอบแบบแม่นตรง (Exact Solution) หรือวิศวกรที่กำลังแก้ปัญหาบางอย่างที่ประกอบไปด้วระบบสมการเชิงเส้นหลายตัวแปร สุดท้ายก็คือให้โปรแกรมเมอร์เห็นถึงวิธีการใช้งานคลาส CMatrix เพื่อนำไปประยุกต์ใช้ในการประมวลผลที่โปรแกรมประยุกต์อื่นๆ มาเริ่มกันเลยครับ

ในตัวอย่างโปรแกรมนี้จะใช้วิธีการโหลดสมาชิกของเมตริกซ์จาก Text file ที่เก็บข้อมูลไว้ ดังที่เคยนำเสนอวิธีการไว้ครับ และเป็นที่ทราบกันดีว่าระบบสมการเชิงเส้นสามารถเขียนในรูปแบบเมตริกซ์ได้ดังนี้

[A][x] = [b]

โดย A คือเมตริกซ์จัตตุรัสขนาด n x n ที่ประกอบด้วยค่าสัมประสิทธิ์คงที่ของตัวแปร x ในแต่ละสมการ

x คือ เมตริกซ์ตัวแปรที่ต้องหาคำตอบขนาด n x 1

b คือ เมตริกซ์ค่าด้านขวาของระบบสมการขนาด n x 1

ตัวอย่างเช่น พิจารณาระบบสมการเส้น 3 ตัวแปร

x + y + z = 10

x - y + z = 6

2x + y - 3z = -17.5

โดยจะบันทึกสมาชิกของเมตริก A และ b ในไฟล์ A.txt และ b.txt ไว้ที่ c:// ตามบทความที่เคยนำเสนอการเขียนโปรแกรมอ่านข้อมูลจาก text ไฟล์ 

มาดู code program ภาษา c++ กันครับ ดูแนวคิดกันก่อนนะครับ

สร้าง object A ของ class CMatrix โดยเรียก constructor สำหรับโหลดข้อมูลจากไฟล์ A.txt

สร้าง object b ของ class CMatrix โดยเรียก constructor สำหรับโหลดข้อมูลจากไฟล์ b.txt

สร้าง object x ของ class CMatrix และเรียกใช้ฟังก์ชัน SolveGuass เพื่อหาคำตอบของระบบสมการ โดยใช้ระเบียบวิธี  Gauss elimination หรือ LU decomposition

แสดงผลการหาคำตอบของระบบสมการเชิงเส้นมาดู Code กันเลยครับ

มาดูการใช้งานโปรแกรมนี้กันเลยครับตามภาพที่ 1

ท่านสามารถโหลด Program LES ไปใช้งานได้เลยนะครับ ตามวัตถุประสงค์ของบทความนี้ แต่ต้องระวังข้อผิดพลาดจากผู้ใช้ เช่น เป็นระบบสมการที่ไม่มีคำตอบ ซึ่งจะทำให้โปรแกรมทำงานผิดพลาด การไม่มีอยู่ของไฟล์เมตริกซ์ A หรือ b ตาม path ที่กำหนด ซึ่งตรงนี้ผมยังไม่ได้ทำ code ตรวจสอบไว้นะครับ ในบทต่อๆไปจะพยายามแทรก code การตรวจสอบไว้ให้เป็นตัวอย่างนะครับ นำเสนอการประยุกต์ใช้ Class CMatrix ในการสร้างโปรแกรมแก้ปัญหาอื่นๆอีกครับ