# linear elliptic vector equation

Barry Smith bsmith at mcs.anl.gov
Wed Aug 19 13:01:10 CDT 2009

```On Aug 19, 2009, at 2:04 AM, nicolas aunai wrote:

> Hi,
>
> OK then, lets go for multigrid.
>
>
> When I'm looking for example at exercice 32 :
> http://www.mcs.anl.gov/petsc/petsc-as/snapshots/petsc-current/src/ksp/ksp/examples/tutorials/ex32.c.html
>
> I realize that I don't understand some things, in particular about the
> DMMGSetKSP function.
>
> I need to give the RHS function and the ComputeMatrix function.
>
> 1/ My operator is always the same during the simulation, I was
> thinking to build the matrix once and for all at the begining of the
> program. Why can't I do that ? Is it because at each multigrid level
> the required matrix is different (size) and it can only be built with
> the information the DA is giving ?

Yes. The idea is rather than you giving the matrix. You give the
function that fills up the matrix.
>
> 2/ I don't really understand the prototype of the callback function
> ComputeMatrix. The first argument is the DMMG object, ok. But the
> second and third arguments are matrices called 'J' and 'jac'.  It
> seems (at least in ex32) that 'J' is not used, what is the goal of
> this matrix in general ?

Sometimes one uses say a second order discretization to define the
system but wants a cheap preconditioner built from the first order
discretization. Hence the two matrices. If you are not doing this then
you pass the same matrix in both slots.

>
> Is 'jac' stands for jacobian ? Reading the ComputeMatrix() code, it
> seems that the 'jac' matrix is built like a standard finite difference
> matrix, with the stencil corresponding to the operator I have to
> solve. if so, why is it called 'jacobian' ? (sorry it looks more a
> math question than a petsc question...)

This is just brought over from the nonlinear solver code. It could
just be called matrix here.

Barry

>
>
> Another quick question :
>
> Why is there no documentation for the function 'DMMGGetr(dmmg)' which
> seems to be part of the petsc library since I can't see the definition
> of the function in the exercice file. I suppose it computes the
> residual ?
>
>
> I guess this is it.... for the moment.
> Thx
>
> Nico
>
>
>
> 2009/8/19 Barry Smith <bsmith at mcs.anl.gov>:
>>
>> On Aug 18, 2009, at 11:42 AM, nicolas aunai wrote:
>>
>>> Hello,
>>>
>>> Yes it is a bit overkill indeed :-) I was looking at ex29
>>> actually...
>>> which is still linear and with DA's but scalar.
>>>
>>> But I'd rather like not use multigrids, however I have not seen an
>>> example of a code that uses DA's and KSP solvers.
>>>
>>> Every example use multigrids.
>>
>>   Multigrid is the way to go for this problem, it will be faster
>> than some
>> generic preconditioned CG.
>>
>>>
>>> When you have a multicomponent linear equation to solve on a
>>> rectangular uniform grid such as mine, is it good to :
>>>
>>> 1/  use DA's to represent your grid
>>>
>>> 2/ create the matrix corresponding to your operator once and for all
>>> (the operator neve changes).( In fact there are 2 matrices because
>>> two
>>> cases possibles (Neumann/periodic or Dirichlet/Periodic, depending
>>> on
>>> the component). but never minds)
>>>
>>>
>>> 3/ Call KSPSolve() 3 times (for each vector component associated
>>> with the
>>> DA) ?
>>>
>>   Yes. You can do all this using the DMMG object.
>>
>>  Note that the DMMG object in your linear case just manages a KSP
>> and the
>> meshes for you. In the end it is a KSPSolve() that actually does
>> the solve
>> and you have full control over the solve including using -pc_type
>> icc or
>> something else that does not use multigrid. Though multigrid will
>> be faster.
>>
>>
>>   Barry
>>
>>
>>>
>>> Thx
>>> Nico
>>>
>>>
>>> 2009/8/18 Aron Ahmadia <aron.ahmadia at kaust.edu.sa>:
>>>>
>>>> ex19 comes to mind, though it's a bit overkill for what you're
>>>> doing...
>>>>
>>>>
>>>> http://www.mcs.anl.gov/petsc/petsc-as/snapshots/petsc-current/src/snes/examples/tutorials/ex19.c.html
>>>>
>>>> ex19 uses DAs and DMMG, which is kind of like a meta-DA for using
>>>> multigrid-style solvers.  Both work well with structured grids.
>>>>
>>>> You may also want to review the PETSc manual section on DA.
>>>>
>>>> A
>>>>
>>>> On Tue, Aug 18, 2009 at 3:30 AM, nicolas aunai <nicolas.aunai at gmail.com
>>>> >
>>>> wrote:
>>>>>
>>>>> Dear all,
>>>>>
>>>>>
>>>>> I have a simulation code (particle in cell) in which I have to
>>>>> solve a
>>>>> linear vector equation 4 times per time step. A typical run
>>>>> consists
>>>>> of 50 000 time steps. The grid is rectangular and uniform of
>>>>> size Lx,
>>>>> Ly, with nx+1 and ny+1 points in the x and y direction
>>>>> respectively.,
>>>>> (nx, ny) could be at max (1024,1024).
>>>>>
>>>>> The vector equation is the following :
>>>>>
>>>>> B(x,y) - a*Laplacian(B(x,y)) = S(x,y)
>>>>>
>>>>>
>>>>> Where B and S are 2D vector fields with 3 components (Bx, By, Bz
>>>>> and
>>>>> Sx, Sy, Sz), each depending on the x and y coordinates.
>>>>>
>>>>>
>>>>> 'a' is a positive constant, smaller than one, typically a= 0.02
>>>>>
>>>>> Boundary conditions may depend on for which B component we are
>>>>> solving
>>>>> the equation. On the x=cst borders, the boundary is always
>>>>> periodic,
>>>>> no matter what component is solved, but on y=cst borders, the
>>>>> boundary
>>>>> condition can be either Neumann or Dirichlet.
>>>>>
>>>>> So far, I have written a small code that creates a vector
>>>>> solution, a
>>>>> vector RHS and the matrix operator, and solve a scalar equation of
>>>>> this type. Solving my vector equation would then just call 3 times
>>>>> this kind of code... but I believe there is another way to deal
>>>>> with
>>>>> vector fields and linear systems with Petsc, using DAs no ? Could
>>>>> someone explain me how solve this the proper way and/or show me
>>>>> some
>>>>> code solving linear system with vector fields ? (is there an
>>>>> example
>>>>> in the exercices that I would not have seen ?)
>>>>>
>>>>>
>>>>>
>>>>> Thanks a lot
>>>>> Nico
>>>>
>>>>
>>>>
>>>> --
>>>> Aron Jamil Ahmadia
>>>> Assistant Research Scientist
>>>> King Abdullah University of Science and Technology
>>>>
>>
>>

```