[petsc-users] Fwd: Re: Problem with fortran version of ex29 in ksp
TAY wee-beng
zonexo at gmail.com
Fri May 4 07:59:00 CDT 2012
Hi,
Is there anything else I can try to get it working right?
Thanks!
On 2/5/2012 10:11 PM, Matthew Knepley wrote:
> On Wed, May 2, 2012 at 1:55 PM, TAY wee-beng <zonexo at gmail.com
> <mailto:zonexo at gmail.com>> wrote:
>
> Hi,
>
> I did a MatView and VecView on both C and Fortran, right after Mat
> and Vec assembly. I have attached the printout below. They are
> exactly the same, but yet the result is different in Neumann
> condition. However, the dirichlet condition gives the correct ans.
> Is there anything else that could be wrong even if the Mat and Vec
> are the same?
>
>
> Did you set the null space for the matrix when you have Neumann
> conditions?
Yes, for the matrix, I set as:
call
MatNullSpaceCreate(PETSC_COMM_WORLD,PETSC_TRUE,0,PETSC_NULL_OBJECT,nullspace,ierr)
call MatSetNullSpace(jac,nullspace,ierr)
call MatNullSpaceDestroy(nullspace,ierr)
for the Vec,
call
MatNullSpaceCreate(PETSC_COMM_WORLD,PETSC_TRUE,0,PETSC_NULL_OBJECT,nullspace,ierr)
*!call MatNullSpaceRemove(nullspace,b,PETSC_NULL,ierr)*
call MatNullSpaceDestroy(nullspace,ierr)
MatNullSpaceRemove was comment out because there's error during linking
>
> Matt
>
> Thanks!
>
> Fortran:
>
> Matrix Object: 1 MPI processes
> type: seqaij
> row 0: (0, 2) (1, -1) (3, -1)
> row 1: (0, -1) (1, 3) (2, -1) (4, -1)
> row 2: (1, -1) (2, 2) (5, -1)
> row 3: (0, -1) (3, 3) (4, -1) (6, -1)
> row 4: (1, -1) (3, -1) (4, 4) (5, -1) (7, -1)
> row 5: (2, -1) (4, -1) (5, 3) (8, -1)
> row 6: (3, -1) (6, 2) (7, -1)
> row 7: (4, -1) (6, -1) (7, 3) (8, -1)
> row 8: (5, -1) (7, -1) (8, 2)
> Vector Object:Vec_0000000084000000_0 1 MPI processes
> type: mpi
> Process [0]
> 0.25
> 0.0205213
> 1.135e-005
> 0.0205213
> 0.00168449
> 9.31663e-007
> 1.135e-005
> 9.31663e-007
> 5.15289e-010
> Vector Object:Vec_0000000084000000_1 1 MPI processes
> type: mpi
> Process [0]
> 0.14924
> 0.0242397
> -0.0260347
> 0.0242397
> -0.0256192
> -0.0400102
> -0.0260347
> -0.0400102
> -0.0400102
> Press any key to continue . . .
>
> C:
>
> Matrix Object: 1 MPI processes
> type: seqaij
> row 0: (0, 2) (1, -1) (3, -1)
> row 1: (0, -1) (1, 3) (2, -1) (4, -1)
> row 2: (1, -1) (2, 2) (5, -1)
> row 3: (0, -1) (3, 3) (4, -1) (6, -1)
> row 4: (1, -1) (3, -1) (4, 4) (5, -1) (7, -1)
> row 5: (2, -1) (4, -1) (5, 3) (8, -1)
> row 6: (3, -1) (6, 2) (7, -1)
> row 7: (4, -1) (6, -1) (7, 3) (8, -1)
> row 8: (5, -1) (7, -1) (8, 2)
> Vector Object:Vec_0x1d3b000_0 1 MPI processes
> type: mpi
> Process [0]
> 0.25
> 0.0205212
> 1.135e-05
> 0.0205212
> 0.00168449
> 9.31663e-07
> 1.135e-05
> 9.31663e-07
> 5.15288e-10
> Vector Object:Vec_0x1d3b000_1 1 MPI processes
> type: mpi
> Process [0]
> 0.139311
> 0.0305751
> -0.0220633
> 0.0305751
> -0.0135158
> -0.042185
> -0.0220633
> -0.042185
> -0.058449
>
>
>
> Yours sincerely,
>
> TAY wee-beng
>
>
> On 1/5/2012 11:54 PM, Matthew Knepley wrote:
>> On Tue, May 1, 2012 at 5:48 PM, TAY wee-beng <zonexo at gmail.com
>> <mailto:zonexo at gmail.com>> wrote:
>>
>> Hi,
>>
>> Do you mean my method is wrong?
>>
>> I am following the template of ex22f,
>>
>> where the variables are declared as :
>>
>> PetscScalar v(5)
>>
>> MatStencil row(4),col(4,5)
>>
>> Hence,
>>
>> for the neumann BC
>>
>> num = 1
>>
>> if (j/=0) then
>>
>> v(num) = -rho*HxdHy
>>
>> col(MatStencil_i,num) = i
>>
>> col(MatStencil_j,num) = j-1
>>
>> num = num + 1
>>
>> end if
>>
>> if (i/=0) then
>>
>> v(num) = -rho*HydHx
>>
>> col(MatStencil_i,num) = i-1
>>
>> col(MatStencil_j,num) = j
>>
>> num = num + 1
>>
>> end if
>>
>> if (i/=mx-1) then
>>
>> v(num) = -rho*HydHx
>>
>> col(MatStencil_i,num) = i+1
>>
>> col(MatStencil_j,num) = j
>>
>> num = num + 1
>>
>> end if
>>
>> if (j/=my-1) then
>>
>> v(num) = -rho*HxdHy
>>
>> col(MatStencil_i,num) = i
>>
>> col(MatStencil_j,num) = j+1
>>
>> num = num + 1
>>
>> end if
>>
>> v(num) = ((num-1)/2.0)*rho*(HxdHy + HydHx)
>>
>> print *, v
>>
>> col(MatStencil_i,num) = i
>>
>> col(MatStencil_j,num) = j
>>
>> !num = num + 1
>>
>> call
>> MatSetValuesStencil(jac,i1,row,num,col,v,INSERT_VALUES,ierr)
>>
>> I do not get any more out of range error. However,my ans is
>> still different from that of ex29 in C.
>>
>>
>> This is very simple. You have an error in your code. Checking it
>> is very simple: run the code and
>> break in MatSetValues(). Make sure ex29 makes calls with exactly
>> the same indices as your ex29f.
>>
>> Matt
>>
>> Yours sincerely,
>>
>> TAY wee-beng
>>
>>
>> --
>> 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
>
>
>
>
> --
> 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
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20120504/51e90ef3/attachment.htm>
More information about the petsc-users
mailing list