[petsc-users] Insert values into matrix using MatSetValuesStencil
TAY wee-beng
zonexo at gmail.com
Fri Sep 25 03:20:42 CDT 2015
Hi Timothée,
I got the error sorted out. I was inserting into the wrong matrix. But
the key error was the defn of row and col.
Thank you for pointing it out.
Yours sincerely,
TAY wee-beng
On 23/9/2015 5:25 PM, Timothée Nicolas wrote:
> Hum, sorry, I don't know. I asked you to provide the definitions of
> start and end values of i,j,k, because I was concerned whether you may
> mess up the boundaries. Especially because you seem to treat x
> differently from y and z. You have the problem also on only 1 process ?
>
> 2015-09-23 18:13 GMT+09:00 TAY wee-beng <zonexo at gmail.com
> <mailto:zonexo at gmail.com>>:
>
> Hi Timothée,
>
> Maybe I can send you part of it 1st. I'm trying to pinpoint why my
> matrix using MatView shows zero for a lot of the values
>
> For i=1,j=1,k=1,
>
> It should be :
>
> Mat Object: 1 MPI processes
> type: seqaij
> row 0: (0, 2) (12, -2)
>
> but now it's:
>
> row 0: (0, 0) (1, 0) (2, 0) (3, 0) (4, 0) (5, 0) (6, 0) (7,
> 0) (8, 0) (12, 0) (13, 0) (14, 0) (24, 0) (25, 0) (26, 0)
> (96, 0) (97, 0) (98, 0) (192, 0) (193, 0) (194, 0)
>
> I used:
>
> /call
> DMDACreate3d(MPI_COMM_WORLD,DM_BOUNDARY_NONE,DM_BOUNDARY_NONE,DM_BOUNDARY_NONE,DMDA_STENCIL_STAR,size_x,size_y,&//
> //
> //size_z,1,PETSC_DECIDE,PETSC_DECIDE,3,stencil_width,lx,PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,da_uvw_ast,ierr)//
> //
> //call DMSetMatType(da_uvw_ast,MATAIJ,ierr)//
> ////
> // call DMCreateMatrix(da_uvw_ast,A_semi_xyz,ierr)//
> //
> // call MatSetFromOptions(A_semi_xyz,ierr)//
> //
> //i = 1, j = 1, k = 1//
> //
> //row(MatStencil_i,1) = i - 1//
> ////
> // row(MatStencil_j,1) = j - 1//
> ////
> // row(MatStencil_k,1) = k - 1//
> ////
> // row(MatStencil_c,1) = 0//
> ////
> // col(MatStencil_i,1) = i - 1//
> ////
> // col(MatStencil_j,1) = j - 1//
> ////
> // col(MatStencil_k,1) = k - 1//
> ////
> // col(MatStencil_c,1) = 0//
> ////
> // value_insert(1) = 2.//
> //
> // call
> MatSetValuesStencil(A_mat,ione,row,ione,col(:,1),value_insert(1),INSERT_VALUES,ierr)//
> //
> // col(MatStencil_i,2) = i - 1 !north//
> ////
> // col(MatStencil_j,2) = j + 1 - 1//
> ////
> // col(MatStencil_k,2) = k - 1//
> ////
> // col(MatStencil_c,2) = 0//
> ////
> // value_insert(2) = -2.//
> //
> // call
> MatSetValuesStencil(A_mat,ione,row,ione,col(:,2),value_insert(2),INSERT_VALUES,ierr)/
>
> Thank you
>
> Yours sincerely,
>
> TAY wee-beng
>
> On 23/9/2015 4:45 PM, Timothée Nicolas wrote:
>> Yes, I had understood that, I am doing the same, but with 8 dof.
>> This does not change the declaration for row and col.
>>
>> Can you send (i) the commands you use to create the DMDA, (ii)
>> the commands to create the matrix and (iii) those for the
>> definitions of ksta2,kend2,jsta2,jend2,size_x ?
>>
>> Timothée
>>
>> 2015-09-23 17:37 GMT+09:00 TAY wee-beng <zonexo at gmail.com
>> <mailto:zonexo at gmail.com>>:
>>
>> Hi Timothée,
>>
>> The matrix is created with 3 dof - u,v,w. So for each i,j,k,
>> there are 3 values. Actually I got 3 eqns from the u,v,w
>> momentum eqns. They are not coupled together so I can also
>> solve them individually. But I was told it's faster to group
>> them together. The global indices is something like this:
>>
>> i, j, k, ijk global indices, dof
>> 1 1 1 1 u
>> 1 1 1 2 v
>> 1 1 1 3 w
>> 2 1 1 4 u
>>
>> ...
>>
>> Hope it's clearer now.
>>
>> Ok, I changed the stencil, it's now working.
>>
>> Thank you
>>
>> Yours sincerely,
>>
>> TAY wee-beng
>>
>> On 23/9/2015 4:22 PM, Timothée Nicolas wrote:
>>> Can you also tell how you created the matrix ? Just in case
>>> you created it with the 1 dof DMDA, it would not work if you
>>> try to input values at places where it is not allocated
>>> (which could explain the error message)
>>>
>>> 2015-09-23 17:18 GMT+09:00 Timothée Nicolas
>>> <timothee.nicolas at gmail.com
>>> <mailto:timothee.nicolas at gmail.com>>:
>>>
>>> The first thing that strikes me is your definition of
>>> the stencils
>>>
>>> */MatStencil :: row(6,1),col(6,7)/**/
>>>
>>> /*
>>> Why is it not defined with
>>>
>>> */MatStencil :: row(4,1),col(4,7)/**/
>>> /*
>>> instead ?
>>>
>>> Where does the 6 come from ?
>>>
>>> Timothée
>>>
>>>
>>> 2015-09-23 17:14 GMT+09:00 TAY wee-beng
>>> <zonexo at gmail.com <mailto:zonexo at gmail.com>>:
>>>
>>> Hi,
>>>
>>> I have successfully used MatSetValuesStencil to
>>> insert values into a Poisson eqn matrix which has 1
>>> DOF (pressure). Now I'm trying to insert values in a
>>> momentum eqn matrix which has 3 DOF (u,v,w)
>>>
>>> However, I got the error:
>>>
>>> /*[0]PETSC ERROR: --------------------- Error
>>> Message ----------------------------*//*
>>> *//*----------------------------------*//*
>>> *//*[0]PETSC ERROR: Argument out of range*//*
>>> *//*[0]PETSC ERROR: Inserting a new nonzero at
>>> (111,5) in the matrix*//*
>>> *//*[0]PETSC ERROR: See
>>> http://www.mcs.anl.gov/petsc/documentation/faq.html
>>> for trou*//*
>>> *//*ble shooting.*/
>>>
>>> I wonder what's wrong. For the momentum eqn, for
>>> each DOF, at at node (dof,i,j,k), I have coupling i
>>> +/- 1, j +/- 1 and k +/- 1.
>>>
>>> The error happens at 111,5, which corresponds to i =
>>> 2, j = 2, k = 2, which is an internal node.
>>>
>>> Here's part of my code below. Hope someone can help.
>>> Thanks!
>>> */
>>> /**/PetscScalar :: value_insert(7)/**/
>>> /**/
>>> /**/MatStencil :: row(6,1),col(6,7)/**/
>>> /**/
>>> /**/ione = 1; iseven = 7/**/
>>> /**/
>>> /**/if (cell_type == 'u') then/**/
>>> /**/
>>> /**/offset = 1/**/
>>> /**//**/
>>> /**/else if (cell_type == 'v') then/**/
>>> /**/
>>> /**/offset = 2/**/
>>> /**//**/
>>> /**/else if (cell_type == 'w') then/**/
>>> /**/
>>> /**/offset = 3/**/
>>> /**//**/
>>> /**/end if/**/
>>> /**/
>>> /**/do k=ksta2,kend2/**/
>>> /**/
>>> /**/do j = jsta2,jend2/**/
>>> /**/
>>> /**/ do i=2,size_x-1/**/
>>> /**//**/
>>> /**/row(MatStencil_i,1) = i - 1/**/
>>> /**//**/
>>> /**/row(MatStencil_j,1) = j - 1/**/
>>> /**//**/
>>> /**/row(MatStencil_k,1) = k - 1/**/
>>> /**//**/
>>> /**/row(MatStencil_c,1) = offset - 1/**/
>>> /**//**/
>>> /**/value_insert = 0.d0/**/
>>> /**//**/
>>> /**/col(MatStencil_i,3) = i + 1 - 1 !east/**/
>>> /**//**/
>>> /**/col(MatStencil_j,3) = j - 1/**/
>>> /**//**/
>>> /**/col(MatStencil_k,3) = k - 1/**/
>>> /**//**/
>>> /**/col(MatStencil_c,3) = offset - 1/**/
>>> /**//**/
>>> /**/value_insert(3) = -(
>>> 1./(cell_x(i)%pd_E+cell_x(i+1)%pd_W))*(c_yz(j,k)%fc_E)*inv_Re/**/
>>> /**//**/
>>> /**/col(MatStencil_i,5) = i - 1 - 1 !west/**/
>>> /**//**/
>>> /**/col(MatStencil_j,5) = j - 1/**/
>>> /**//**/
>>> /**/col(MatStencil_k,5) = k - 1/**/
>>> /**//**/
>>> /**/col(MatStencil_c,5) = offset - 1/**/
>>> /**//**/
>>> /**/value_insert(5) = -(
>>> 1./(cell_x(i)%pd_W+cell_x(i-1)%pd_E))*(c_yz(j,k)%fc_E)*inv_Re/**/
>>> /**//**/
>>> /**/col(MatStencil_i,2) = i - 1 !north/**/
>>> /**//**/
>>> /**/col(MatStencil_j,2) = j + 1 - 1/**/
>>> /**//**/
>>> /**/col(MatStencil_k,2) = k - 1/**/
>>> /**//**/
>>> /**/col(MatStencil_c,2) = offset - 1/**/
>>> /**//**/
>>> /**/value_insert(2) = -(
>>> 1./(cell_y(j)%pd_N+cell_y(j+1)%pd_S))*(c_zx(i,k)%fc_N)*inv_Re/**/
>>> /**//**/
>>> /**/col(MatStencil_i,4) = i - 1 !south/**/
>>> /**//**/
>>> /**/col(MatStencil_j,4) = j - 1 - 1/**/
>>> /**//**/
>>> /**/col(MatStencil_k,4) = k - 1/**/
>>> /**//**/
>>> /**/col(MatStencil_c,4) = offset - 1/**/
>>> /**//**/
>>> /**/value_insert(4) = -(
>>> 1./(cell_y(j)%pd_S+cell_y(j-1)%pd_N))*(c_zx(i,k)%fc_N)*inv_Re/**/
>>> /**//**/
>>> /**/col(MatStencil_i,6) = i - 1 !front/**/
>>> /**//**/
>>> /**/col(MatStencil_j,6) = j - 1/**/
>>> /**//**/
>>> /**/col(MatStencil_k,6) = k + 1 - 1/**/
>>> /**//**/
>>> /**/col(MatStencil_c,6) = offset - 1/**/
>>> /**//**/
>>> /**/value_insert(6) = -(
>>> 1./(cell_z(k)%pd_F+cell_z(k+1)%pd_B))*(c_xy(i,j)%fc_F)*inv_Re/**/
>>> /**//**/
>>> /**/col(MatStencil_i,7) = i - 1 !back/**/
>>> /**//**/
>>> /**/col(MatStencil_j,7) = j - 1/**/
>>> /**//**/
>>> /**/col(MatStencil_k,7) = k - 1 - 1/**/
>>> /**//**/
>>> /**/col(MatStencil_c,7) = offset - 1/**/
>>> /**//**/
>>> /**/value_insert(7) = -(
>>> 1./(cell_z(k)%pd_B+cell_z(k-1)%pd_F))*(c_xy(i,j)%fc_F)*inv_Re/**/
>>> /**//**/
>>> /**/col(MatStencil_i,1) = i - 1/**/
>>> /**//**/
>>> /**/col(MatStencil_j,1) = j - 1/**/
>>> /**//**/
>>> /**/col(MatStencil_k,1) = k - 1/**/
>>> /**//**/
>>> /**/col(MatStencil_c,1) = offset - 1/**/
>>> /**//**/
>>> /**/value_insert(1) = 2.*c(i,j,k)%vol/del_t -
>>> (value_insert(2)+value_insert(3)+value_insert(4)+value_insert(5)+value_insert(6)+value_insert(7))/**/
>>> /**//**/
>>> /**/call
>>> MatSetValuesStencil(A_semi_xyz,ione,row,iseven,col,value_insert,INSERT_VALUES,ierr)/**/
>>> /**//**/
>>> /**/end do/**/
>>> /**//**/
>>> /**/end do/**/
>>> /**/
>>> /**/end do /*
>>>
>>> Thank you
>>>
>>> Yours sincerely,
>>>
>>> TAY wee-beng
>>>
>>>
>>>
>>>
>>
>>
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20150925/277267ca/attachment-0001.html>
More information about the petsc-users
mailing list