# linear elliptic vector equation

nicolas aunai nicolas.aunai at gmail.com
Wed Aug 19 02:04:19 CDT 2009

```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 ?

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 ?

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...)

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
>>
>>
>>>
>>> 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
>>>
>>>
>>>
>>> --