Here you have an example of a 2 dof system using MatSetValuesBlockedStencil
https://github.com/stefanozampini/petscopt/blob/master/src/ts/examples/tests/ex3.c#L186 <https://github.com/stefanozampini/petscopt/blob/master/src/ts/examples/tests/ex3.c#L186> > On May 29, 2020, at 8:55 PM, Lucas Banting <banti...@myumanitoba.ca> wrote: > > Yes, I think that makes sense to me. Each cell is at i,j but then you can set > the 25 coefficients per node using MatStencil_c from 0 to 4 for both idxm and > idxn. > > Thanks for your help. > > Lucas > From: Matthew Knepley <knep...@gmail.com> > Sent: Friday, May 29, 2020 12:47 PM > To: Lucas Banting <banti...@myumanitoba.ca> > Cc: petsc-users@mcs.anl.gov <petsc-users@mcs.anl.gov> > Subject: Re: [petsc-users] How to use MatSetValuesStencil() with DOF>1 > > On Fri, May 29, 2020 at 12:16 PM Lucas Banting <banti...@myumanitoba.ca > <mailto:banti...@myumanitoba.ca>> wrote: > Hello, > > I have a structured grid problem with 5 unknowns per cell: U,V,P,T,F, this is > a CFD code but I don't think these details are necessary to answer my > question. > I was wondering how exactly I should use MatSetValuesStencil() to fill my > matrix. > I have 5 degrees of freedom, so I need to set up to 25 coefficients per grid > cell. I don't understand how i,j, and c map to all coefficients. > > For a 1 degree of freedom system to fill my 9 point stencil I did: > > values(1)= ASW(II,JJ) > values(2)= AS (II,JJ) > values(3)= ASE(II,JJ) > values(4)= AW (II,JJ) > values(5)= AP (II,JJ) > values(6)= AE (II,JJ) > values(7)= ANW(II,JJ) > values(8)= AN (II,JJ) > values(9)= ANE(II,JJ) > > idxm(MatStencil_i,1) = II;idxm(MatStencil_j,1) = JJ > > idxn(MatStencil_i,1) = II-1;idxn(MatStencil_j,1) = JJ-1 > idxn(MatStencil_i,2) = II ;idxn(MatStencil_j,2) = JJ-1 > idxn(MatStencil_i,3) = II+1;idxn(MatStencil_j,3) = JJ-1 > idxn(MatStencil_i,4) = II-1 ;idxn(MatStencil_j,4) = JJ > idxn(MatStencil_i,5) = II ;idxn(MatStencil_j,5) = JJ > idxn(MatStencil_i,6) = II+1;idxn(MatStencil_j,6) = JJ > idxn(MatStencil_i,7) = II-1 ;idxn(MatStencil_j,7) = JJ+1 > idxn(MatStencil_i,8) = II ;idxn(MatStencil_j,8) = JJ+1 > idxn(MatStencil_i,9) = II+1;idxn(MatStencil_j,9) = JJ+1 > call MatSetValuesStencil(A,1,idxm,9,idxn,values,INSERT_VALUES,ierr) > > Which seemed to work just fine. For 5 degrees of freedom, instead of just > having a single AP coefficient for example, I have an AP matrix: > > ap_uu > ap_uv > ap_up > ap_ut > ap_uf > ap_vu > ap_vv > ap_vp > ap_vt > ap_vf > ap_pu > ap_pv > ap_pp > ap_pt > ap_pf > ap_tu > ap_tv > ap_tp > ap_tt > ap_tf > ap_fu > ap_fv > ap_fp > ap_ft > ap_ff > > In 1 degree of freedom, i and j corresponded to the solution variable in the > grid. For multi degree of freedom, I don't understand how values of i, j, and > c could distinguish from an ap_uv and ap_vv coefficent for example, wouldn't > they both be at i, j ,c=2? > Is there anyway I can use MatSetValuesStencil() to fill in my 9 point stencil > with my 5x5 matrix coefficients? To clarify I have nine 5x5 matrices for each > cell which correspond to the 5 unknowns per cell. > > In MatSetValuesStencil(), I think it is best to think of the Stencils as > another way of providing row/col numbers. So each MatStencil corresponds > to some row number. So your ap_uv is some entry in the input matrix, and thus > corresponds to some row number (a MatStencil) and some > column number (another MatStencil). Since ap is one point, it appears that > > ap_uv --> row = i,j,0 col = i,j,1 > ap_vv --> row = i,j,1 col = i,j,1 > > So you line up your MatStencil arguments to match the order of your input. It > sounds like you want the field index to be the fastest in your input. > > Does that make sense? > > Thanks, > > Matt > > > Thanks, > > Lucas > > > -- > 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 > > https://www.cse.buffalo.edu/~knepley/ <http://www.cse.buffalo.edu/~knepley/>