ONLY DO WHAT ONLY YOU CAN DO

こけたら立ちなはれ 立ったら歩きなはれ

Go言語 で連立一次方程式を解く(ヤコビの反復法)

例として
9x+2y+z+u=20
2x+8y-2z+u=16
-x-2y+7z-2u=8
x-y-2z+6u=17
を考える.
この方程式を上から順に対角線上の変数について解くと
x=(20-2y-z-u)/9
y=(16-2x+2z-u)/8
z=(8+x+2y+2u)/7
u=(17-x+y+2z)/6
となる.
x,y,z,uに適当な値を入れて右辺を計算し、
得られた値を新たなx,y,z,uとして、計算を繰り返す.
漸化式で書くと
x_{n+1}=(20-2y_n-z_n-u_n)/9
y_{n+1}=(16-2x_n+2z_n-u_n)/8
z_{n+1}=(8+x_n+2y_n+2u_n)/7
u_{n+1}=(17-x_n+y_n+2z_n)/6

package main

import "fmt"
import "math"

const N = 4

func main() {
    var a [N][N]float64 = [N][N]float64{{9,2,1,1},{2,8,-2,1},{-1,-2,7,-2},{1,-1,-2,6}} 
    var b=[]float64{20,16,8,17}
    var c=[]float64{0,0,0,0}

    jacobi(a,b,c)

    fmt.Println("X")
    for i := 0; i < N; i++ {
        fmt.Printf("%14.10f\t", c[i])
    }
    fmt.Println("")
}
func jacobi(a[N][N]float64, b[]float64, c[]float64) {
    for {
        var s=[]float64{0,0,0,0}
        var finish = true
        for i := 0; i < N; i++ {
            s[i] = b[i]
            for j := 0; j < N; j++ {
                if (j != i) {
                    s[i] -= a[i][j] * c[j]
                }
            }
            s[i] /= a[i][i]
            if (math.Abs(s[i] - c[i]) > 0.0000000001) {
              finish = false
            }
        }
        for i := 0; i < N; i++ {
            c[i] = s[i]
            fmt.Printf("%14.10f\t", c[i])
        }
        fmt.Println("")

        if (finish) {
          return
        }
    }
}
  2.2222222222	  2.0000000000	  1.1428571429	  2.8333333333	
  1.3359788360	  1.3759920635	  2.8412698413	  3.1772486772	
  1.2477219283	  1.9791666667	  2.6346371882	  3.7870921517	
  1.0688819252	  1.8733422960	  2.9686056521	  3.8334531858	
  1.0501396189	  1.9957492835	  2.9260675555	  3.9569452792	
  1.0139431776	  1.9743638243	  2.9936469635	  3.9662907960	
  1.0101482880	  1.9991395970	  2.9850360597	  3.9912857623	
  1.0028221093	  1.9948112226	  2.9987141438	  3.9931772381	
  1.0020540192	  1.9998258539	  2.9969712901	  3.9982362335	
  1.0005711965	  1.9989497885	  2.9997397420	  3.9986190691	
  1.0004157346	  1.9999647527	  2.9993869874	  3.9996430127	
  1.0001156105	  1.9997874366	  2.9999473236	  3.9997204988	
  1.0000841449	  1.9999928659	  2.9998759259	  3.9999277456	
  1.0000233996	  1.9999569771	  2.9999893383	  3.9999434288	
  1.0000170310	  1.9999985561	  2.9999748873	  3.9999853757	
  1.0000047361	  1.9999912921	  2.9999978421	  3.9999885500	
  1.0000034471	  1.9999997077	  2.9999949172	  3.9999970400	
  1.0000009586	  1.9999982375	  2.9999995632	  3.9999976825	
  1.0000006977	  1.9999999408	  2.9999989712	  3.9999994009	
  1.0000001940	  1.9999996433	  2.9999999116	  3.9999995309	
  1.0000001412	  1.9999999880	  2.9999997918	  3.9999998787	
  1.0000000393	  1.9999999278	  2.9999999821	  3.9999999051	
  1.0000000286	  1.9999999976	  2.9999999579	  3.9999999755	
  1.0000000079	  1.9999999854	  2.9999999964	  3.9999999808	
  1.0000000058	  1.9999999995	  2.9999999915	  3.9999999950	
  1.0000000016	  1.9999999970	  2.9999999993	  3.9999999961	
  1.0000000012	  1.9999999999	  2.9999999983	  3.9999999990	
  1.0000000003	  1.9999999994	  2.9999999999	  3.9999999992	
  1.0000000002	  2.0000000000	  2.9999999997	  3.9999999998	
  1.0000000001	  1.9999999999	  3.0000000000	  3.9999999998	
  1.0000000000	  2.0000000000	  2.9999999999	  4.0000000000	
  1.0000000000	  2.0000000000	  3.0000000000	  4.0000000000	
X
  1.0000000000	  2.0000000000	  3.0000000000	  4.0000000000	
参考文献