[petsc-users] Solving Ill Conditioned Problems

[Barry Smith]

  I would do the scaling in the model before ever writing any code.  The problem is that when the numbers are of very different sizes you have to be careful with ANY computations (though, as note, adding and subtracting is the most dangerous).

    Barry


On Oct 19, 2012, at 2:09 PM, Nachiket Gokhale <gokhalen at gmail.com> wrote:

> Yes, I am doing the scaling in the element routines for the matrix assembly.
> 
> I had not realized the importance of scaling. I was under the
> impression that if vastly different scales are in different blocks of
> the matrix, then the preconditioners or balancing operations will take
> care of the scaling.  To be more precise, I was under the impression
> that things should be ok as long as one is not adding or subtracting
> numbers of vastly different scales to each other.
> 
> -Nachiket
> 
> On Fri, Oct 19, 2012 at 3:04 PM, Barry Smith <bsmith at mcs.anl.gov> wrote:
>> 
>> On Oct 19, 2012, at 11:02 AM, Anton Popov <popov at uni-mainz.de> wrote:
>> 
>>> Of, course you have a more sophisticated case, but if you just put the orders of magnitude of your matrix blocks
>>> into a toy two-by-two matrix, i.e. (in matlab notation)
>>> A = [1e9, 1; 1, 1e-9]
>>> then it's immediately visible that A is singular (det(A) is zero)
>>> Symmetric scaling, that Matt suggested, just makes it more obvious, after scaling the matrix A will be:
>>> A = [1, 1; 1, 1] - clearly the rows are linearly dependent.
>>> Scaling may help you, but of course, this is not an always working solution
>> 
>>   You need to do the scaling BEFORE you construct the matrix (or ideally before you do computing). Yes, forming the matrix and then scaling it is hopeless
>> 
>>   Barry
>> 
>>> 
>>> Anton
>>> 
>>> On 10/18/12 8:38 PM, Nachiket Gokhale wrote:
>>>> On Thu, Oct 18, 2012 at 2:28 PM, Matthew Knepley <knepley at gmail.com> wrote:
>>>>> On Thu, Oct 18, 2012 at 11:07 AM, Nachiket Gokhale <gokhalen at gmail.com>
>>>>> wrote:
>>>>>> Hi,
>>>>>> 
>>>>>> I am solving a Piezoelectric problem using Petsc. The structure is
>>>>>> 
>>>>>> [ K_uu        K_uv ]
>>>>>> [K_uv ^T   -K_v,v ]
>>>>>> 
>>>>>> More details about the formulation:  http://tinyurl.com/9hlbp4u
>>>>>> 
>>>>>> K_uu  has elements O(1E9)  because the stiffnesses are in GPa,
>>>>>> K_uv   has elements O(1) because piezoelectric coefficients are of that
>>>>>> order
>>>>>> K_v,v  has  elements O(1E-9) because the dielectric constants are of
>>>>>> that order.
>>>>>> 
>>>>>> I am using Petsc, with pc_type LU and MUMPS for the factorization,
>>>>>> -ksp_type gmres. I am not sure if my solution is converging. A typical
>>>>>> solve seems to be doing this:
>>>>>> 
>>>>>> 28 KSP preconditioned resid norm 5.642364260456e-06 true resid norm
>>>>>> 1.228976487745e-03 ||r(i)||/||b|| 3.317409023627e-14
>>>>>> 29 KSP preconditioned resid norm 5.540718271043e-06 true resid norm
>>>>>> 1.228453548651e-03 ||r(i)||/||b|| 3.315997440178e-14
>>>>>> 30 KSP preconditioned resid norm 1.973052106578e-03 true resid norm
>>>>>> 1.220399172500e-03 ||r(i)||/||b|| 3.294256047735e-14
>>>>>> 31 KSP preconditioned resid norm 1.155762663956e-17 true resid norm
>>>>>> 2.447631111938e-04 ||r(i)||/||b|| 6.606955965570e-15
>>>>>> 
>>>>>> Is there a right way to solve this set of equations? Is PCFieldSplit
>>>>>> the recommended way?
>>>>> 
>>>>> First, you should really non-dimensionalize. You can see what this would
>>>>> give
>>>>> you by symmetrically scaling your problem with [ 3e4 3e-4 ], namely
>>>>> everything
>>>>> will be O(1).
>>>>> 
>>>>> Second, you might get something out of using FieldSplit, but its tough to
>>>>> tell
>>>>> without knowing more about the operators.
>>>>> 
>>>>>   Matt
>>>>> 
>>>>>> Thanks,
>>>>>> 
>>>>>> -Nachiket
>>>> Oh, right, thanks. I wasn't even thinking of it that way. I'll scale
>>>> the variables and I'll give it a try.
>>>> 
>>>> Cheers,
>>>> 
>>>> -Nachiket
>>>> 
>>>> 
>>>> 
>>>>> 
>>>>> 
>>>>> --
>>>>> What most experimenters take for granted before they begin their experiments
>>>>> is infinitely more interesting than any results to which their experiments
>>>>> lead.
>>>>> -- Norbert Wiener
>>> 
>>