逆行列を求める関数です。本当は正規方程式は掃き出し法で求めると効率がいいのですが、逆行列を使った方が若干分りやすいかと思います。

# 逆行列を求める関数 inv.awk
function inv(r, ip, epsinv,     i, j, k, b, t, l, m)
{
   for (k = 0; k < ip; k++) {
       m[k] = l[k] = k; b = r[k, k]
       for (j = k; j < ip; j++) {
           for (i = k; i < ip; i++) {
               if (abs(r[j, i]) > abs(b)) {
                   b = r[j, i]
                   l[k] = i
                   m[k] = j
               }
           }
       }
       if (abs(b) < epsinv) {
           return 1
       }
       j = l[k]
       if (j > k) {
           for (i = 0; i < ip; i++) {
               t = -r[i, k]
               r[i, k] = r[i, j]
               r[i, j] = t
           }
       }
       i = m[k]
       if (i > k) {
           for (j = 0; j < ip; j++) {
               t = -r[k, j]
               r[k, j] = r[i, j]
               r[i, j] = t
           }
       }
       for (i = 0; i < ip; i++) {
           if (i != k) {
               r[k, i] /= -b
           }
       }
       for (i = 0; i < ip; i++) {
           for (j = 0; j < ip; j++) {
               if (i != k && j != k) {
                   r[j, i] += r[k, i]*r[j, k]
               }
           }
       }
       for (j = 0; j < ip; j++) {
           if (j != k) {
               r[j, k] /= b
           }
       }
       r[k, k] = 1.0/b
   }
   for (k = ip-2; k >= 0; k--) {
       i = l[k]
       if (i > k) {
           for (j = 0; j < ip; j++) {
               t = r[k, j]
               r[k, j] = -r[i, j]
               r[i, j] = t
           }
       }
       j = m[k]
       if (j > k) {
           for (i = 0; i < ip; i++) {
               t = r[i, k]
               r[i, k] = -r[i, j]
               r[i, j] = t
           }
       }
   }
   return 0
}