[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
inverse Matrix problem, Newton's method
From: |
mecht |
Subject: |
inverse Matrix problem, Newton's method |
Date: |
Sat, 10 Nov 2018 14:35:59 -0600 (CST) |
I'm writing some Octave code to solve a finite element-like problem using the
Newton's method.
My Jacobian matrix looks like this (spy(J)):
<img src="https://preview.ibb.co/hknY5A/spy.jpg"; alt="spy" border="0" />
<https://ibb.co/fCC0kA>
Diagonal elements are negative sums of other elements in a row.
My problem is that i can not produce any result due to false(?) matrix
singularity. I tried left-hand matrix division, LU method, inv(A) but Octave
warns about matrix singularity due to machine precision.
Octave returns 0 or extremely low value (like 1.0E-100) when asked for
det(A)
At first glance, i thought the reason is that if I multiply large amount
diagonal values of 1.0E-1 order the product falls under machine precision.
BUT
I tried my code with minimal number of finite elements 2x2x2, which gives
Jacobian on 8x8 size, and the result was not that low.
Result example:
11.20161 -12.53016 0.00000 0.00000 1.32856 0.00000 0.00000
0.00000
-12.53016 8.54449 0.00000 0.00000 0.00000 3.98567
0.00000 0.00000
0.00000 0.00000 -13.85872 12.53016 0.00000 0.00000
1.32856 0.00000
0.00000 0.00000 12.53016 -16.51584 0.00000 0.00000
0.00000 3.98567
1.32856 0.00000 0.00000 0.00000 -13.85872 12.53016
0.00000 0.00000
0.00000 3.98567 0.00000 0.00000 12.53016 -16.51584
0.00000 0.00000
0.00000 0.00000 1.32856 0.00000 0.00000 0.00000
-13.85872 12.53016
0.00000 0.00000 0.00000 3.98567 0.00000 0.00000
12.53016 -16.51584
Octave says, that det(J):
ans = 0
circshifted diagonals (for positive term of Sarrus' scheme, negative term
looks to be 0):
11.2016 8.5445 -13.8587 -16.5158 -13.8587 -16.5158 -13.8587
-16.5158
0.00000 -12.53016 0.00000 12.53016 0.00000 12.53016
0.00000 12.53016
0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0
1.3286 3.9857 1.3286 3.9857 1.3286 3.9857 1.3286 3.9857
0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0
-12.53016 0.00000 12.53016 0.00000 12.53016 0.00000
12.53016 0.00000
And prod(x) of each diagonal:
1.1477e+09 -0.0000e+00 0.0000e+00 0.0000e+00 7.8619e+02
0.0000e+00 0.0000e+00 -0.0000e+00
>>
Obviously my determinant in this case should be on the order of 1.0E+9
then. Why does this happens?
I'm not experienced much in FEM or CFD code in Octave, but i really wanted
to recreate some simulation from certain publication. How to solve huge set
of equations then? Maybe i'm missing something on maths level here? Do
people even inverse matrices in their finite elements codes then? It is easy
to get
+-Inf or 0 diagonal product if your matrix is large enough. Can (or should)
I increase the precision of computations somehow?
--
Sent from: http://octave.1599824.n4.nabble.com/Octave-General-f1599825.html
- inverse Matrix problem, Newton's method,
mecht <=