אנליזה נומרית/פתרון מערכת משוואות לינאריות

בעיות לינאריות בדרך כלל קלות יותר לפתרון מאשר בעיות לא לינאריות. ניתן לגשת לבעיות אלו בצורה אנליטית, ואז הפתרון תלוי בדיוק המחשב, או בצורה נומרית, ואז הפתרון תלוי במידת התכנסות השיטה וגם בדיוק המחשב. השיטות הנומריות בדרך כלל מהירות יותר מהשיטות האנליטיות (אחרת לא היו משתמשים בהן).

מאחר ואנו מתעניינים בבעיות בעלות פתרון (יחיד), נתקל במטריצות ריבועיות בלבד, אשר מקיימות . לשם נוחות, נציג את אופני הכתיבה השונים של מערכת משוואות לינאריות:

שיטות אנליטיות

עריכה

בשיטות האנליטיות מבצעים פעולות שורה אלמנטריות על המטריצה לקבל הפתרון הרצוי, אשר תלוי בדיוק המחשב עליו מתבצע החישוב. נקודת התורפה המרכזית של השיטות האנליטיות היא גודל המערכת המירבי בו ניתן לטפל.

שיטת גאוס

עריכה

בשיטת גאוס (הנקראת גם "שיטת הדרוג" או "שיטת החילוץ") מדרגים את המטריצה המורחבת [A|b] ובסוף מציבים לאחור. עבור מטריצה מסדר n יידרשו   פעולות.

שיטת גאוס-ג'ורדן

עריכה

עבור מטריצה מסדר n יידרשו   פעולות.

פירוק LU

עריכה

מאלגברה לינארית ידוע, כי כל מטרציה שניתן להביא אותה לצורה משולשת-עליונה על ידי פעולות שורה אלמנטריות, ניתן ליצג באמצעות מכפלת מטריצה תחתונה (L) במטריצה עליונה (U). נדגים זאת באמצעות מטריצה מסדר 3:

 

כעת ניתן למצוא את המקדמים lij, uij כתלות ב-aij. אחרי שמצויות בידינו שתי המטריצות L,U ניתן לפתור בקלות יחסית מערכת משוואות מהצורה Ax=B:

 

מאחר ו-L משושלת תחתונה, ניתן למצוא את y ללא קושי, ואז נשאר לפתור Ux=y, מה שגם מתאפשר בקלות מכיוון ש-U משולשת עליונה.

היפוך מטריצה

עריכה

על מנת לפתור מערכת משוואות מהצורה Ax=B נוכל למצוא את המטריצה ההופכית של A ולבצע מכפלת מטריצות: x=A-1B. בדרך כלל לא נשתמש בשיטה זו כי הביטוי האנליטי למטריצה הופכית מערב מציאת קופקטורים (זהו למעשה ה-Adjoint). נציג את הביטוי בכל זאת:

 

כאשר Cij הם הקופקטורים של המטריצה A, אשר מתקבלים כזכור, על ידי חישוב מינורי המטריצה A.

שיטה אחרת:
בהינתן מערכת מסדר n, נפתור את n הבעיות  , כאשר   הם וקטורי הבסיס הסטנדרטי. לשם פתרון n הבעיות הללו ניתן להעזר בשיטות המופיעות בדף זה.

פתרון מערכת תלת-אלכסונית

עריכה

פתרון מערכת תלת-אלכסונית מתבצע באמצעות אלגוריתם הנקרא TDMA (Tridiagonal matrix algorithm), או אלגוריתם תומס.

 

יש לנו אם כן, 3N-2 איברים, אשר חלקם יכולים להיות 0.

על מנת להגיע לפתרון, נבצע דירוג של שורה אחר שורה: נכפיל את השורה הראשונה ב-   ונוסיף אותה לשורה השנייה, כך שנקבל:

 

כעת נכפיל את את המשוואה ב-   ונוסיף אותה לשורה השלישית, כך שנקבל:

 

כעת, אחרי שעברנו לכתיב β,δ ניתן להכליל ולומר שהמשוואות ה-(i-1)-ית וה-i-ית הן מן הצורה:

 
 

בהתאמה, כך שכאשר נכפיל את המשוואה ה-(i-1)-ית ב-  , ונוסיף למשוואה ה-i-ית, נקבל:

 

לשם השלמת התמונה, נביט בשתי המשוואות האחרונות:

 
 

נכפיל את המשוואה ה-(N-1)-ית ב-  , ונוסיף למשוואה ה-N-ית, כך שנקבל:

 

מכאן מחלצים את xN ומקבלים את שאר הנעלמים על ידי הצבה לאחור.

כאשר נצטרך לכתוב תכנית מחשב, נרצה להכליל את כתיב β,δ גם על השורה הראשונה. לשם כך נוכל להגדיר:

  1.  , ואז האלגוריתם ימשיך:
  2.   עבור  .
  3.  
  4.   עבור  .

בהערכה פשוטה מתקבל כי מספר הפעולות המקסימלי לשיטה זו הינו 5N-4.

שימושים:

  • בהינתן מד"ר, אם נכתוב את משוואת הפרשים עבורה, נקבל מטריצה תלת-אלכסונית.

קישורים חיצוניים

עריכה

שיטות איטרטיביות

עריכה

מאחר והתרגלנו לסמן את מספר האיטרציה ב"n", נסמן את סדר המערכת ב-N.

שיטת Jacobi

עריכה

שיטה זו משתמשת באיברי המטריצה A על מנת להתכנס לפתרון בדרך המהירה ביותר יש לבחור ניחוש התחלתי השווה לממוצע של הערכים העצמיים. מבצעים איטרציות עד להתכנסות של כל וקטור הנעלמים. כלומר: לא מתבצעות איטרציות עבור כל אחד ואחד מהנעלמים בנפרד. מתוך הסכום :  נבודד את הנעלם xi:

 

ואז השיטה האיטרטיבית היא:

 

מה שמתרחש בפועל הוא 3 לולאות מקוננות אשר משתמשות בוקטור   על מנת לייצר את הוקטור   (ראו "קישורים חיצוניים" עבור האלגוריתם).

קריטריוני התכנסות
ניתן להשתמש באחד מן הקריטריונים הבאים:

  •  
  •  
  •  

בדרך כלל התכנסות השיטה היא איטית ולכן יש לבצע מספר רב של איטרציות. את בדיקת ההתכנסות נהוג לבצע בתום הלולאה עבור כל נעלם בנפרד, ולא עבור כל איטרציה בנפרד.

התנאי להתכנסות
אם   הוא הפתרון, אז  . נציב את וקטור השגיאה   לתוך שיטת Jacobi:

 
 

לשם נוחות, נגדיר את השגיאה המקסימלית באיטרציה:   ואז:

 

ואז התנאי להתכנסות הוא:

 

כלומר איברי האלכסון בכל שורה במטריצה A צריכים להיות גדולים מסכום כל שאר האיברים באותה השורה, ואז ההתכנסות מובטחת. ניתן להוכיח שנתאי זה מספיק אך לא הכרחי. לתנאי זה קוראים גם בשם "שליטה אלכסונית".

כעת כשאנו מודעים לתנאי ההתכנסות, ננסה לסדר את שורות המטריצה כך שהתנאי יתקיים, לפני הפעלת השיטה.

שיטת Gauss-Seidel

עריכה

שיטת GS מפצלת את הסכימה לנעלמים לפני הנעלם הנוכחי ולנעלמים אחרי הנעלם הנוכחי:

 

בשיטה זו, בכל איטרציה משתמשים בערכים האחרונים שהתקבלו. כלומר כאן i-1 הנעלמים בוקטור הנעלמים מתעדכנים לפי הסדר יחד עם הנעלם ה-i, עם התקדמות הלולאה. מסיבה זו ההתכנסות מהירה פי 2 משיטת Jacobi. בדיקת ההתכנסות תתבצע כמו בשיטה הקודמת.

השוואה בין שיטת יעקובי לשיטת גאוס-זיידל

עריכה

ננתח את השיטות במקרה של מערכת מסדר 2:

 

במקרה זה, התהליך האיטרטיבי הוא מהצורה:

 

כלומר ההבדל היחיד הוא ששיטת גאוס-זיידל משתמשת בערך המעודכן של x1.

נציב את הביטוי לשגיאה   ונקבל:

 

נזכור כי הוקטור   פותר את המערכת, ואז באמצעות המשוואות הנ"ל נוכל למצוא קשר בין שתי שגיאות עוקבות:

 

כך שעבור שיטת יעקובי מתקיים:

 

ואילו עבור שיטת גאוס-זיידל מתקיים:

 

כלומר בעוד שבשיטת יעקובי σ היא השגיאה כעבור כל שתי איטרציות, בשיטת גאוס-זיידל σ היא השגיאה בכל איטרציה בודדת. לכן בשיטת גאוס-זיידל ההתכנסות מהירה פי 2.

שיטת Successive Over-Relaxation

עריכה

שיטת SOR נועדה על מנת לזרז את התכנסותה של שיטת Gauss-Seidel.

 

ω הוא פרמטר הרלקסציה אשר מקיים  . על מנת להאיץ תהליך התכנסות איטי לוקחים ω>1 ואילו על מנת להבטיח התכנסות עבור שיטות מתבדרות, לוקחים ω<1. שימו לב כי עבור ω=1 מקבלים חזרה את שיטת GS.

כמו כן, עבור מטריצה A נתונה, קיים ערך אופטימלי של ω אשר מביא להתכנסות המהירה ביותר.

ניתן להוכיח שאם שיטת GS מתכנסת, אז שיטת SOR תתכנס מהר יותר.

דרך אחרת לפיתוח השיטה הוא באמצעות שימוש ייצוג המטריצה A באמצעות המכפלה DLU, כאשר D היא מטריצה אלכסונית (למידע נוסף ראו "קישורים חיצוניים").

קישורים חיצוניים

עריכה


הפרק הקודם:
יציבות נומרית
פתרון מערכת משוואות לינאריות הפרק הבא:
פתרון משוואות דיפרנציאליות