[petsc-users] Invite to Fluid Dynamics Software Infrastructure Workshop: April 12-13 at CU Boulder

2019-03-12 Thread Jed Brown via petsc-users
For PETSc users/developers interested in software and data
infrastructure for fluid dynamics:

We are excited to invite you to attend the second workshop to aid in the
conceptualization of FDSI, a potential NSF-sponsored Institute dedicated
to Fluid Dynamics Software Infrastructure. The workshop will be in
Boulder, CO and will start at 1:00 PM on Friday April 12th and continue
until 5:00 PM on April 13th.  Registration and further info about FDSI
is available in the links below.
 
Registration (free; travel support available)

  https://www.colorado.edu/events/cfdsi/fdsi-full-community-workshop
 
GitHub repo of May 2018 Kickoff Workshop:

  https://github.com/FDSI/Kickoff_Workshop
 
Community Needs Assessment based on 2018 Kickoff Workshop:

  https://www.colorado.edu/events/cfdsi/2018-kick-work/Summary
 
FDSI Discussion Forum:

  https://gitter.im/FDSI/community#
 
We also ask for your help in growing the community by either forwarding
this email to any of your colleagues with an interest in fluid dynamics
software or by nominating them for personal invites to this and future
FDSI workshops:

  https://www.colorado.edu/events/cfdsi/content/grow-fdsi
 
We hope to see you in Boulder and/or in online virtual discussions of
fluid dynamics software infrastructure.
 
Thanks,
 
FDSI Conceptualization Management Team


Re: [petsc-users] PetscScatterCreate type mismatch after update.

2019-03-12 Thread Jed Brown via petsc-users
Manuel Valera  writes:

> Ok i'll try that and let you know, for the time being i reverted to 3.9 to
> finish a paper, will update after that :)

3.10 will also work.


Re: [petsc-users] PetscScatterCreate type mismatch after update.

2019-03-12 Thread Manuel Valera via petsc-users
Ok i'll try that and let you know, for the time being i reverted to 3.9 to
finish a paper, will update after that :)

On Tue, Mar 12, 2019 at 6:42 PM Zhang, Junchao  wrote:

> Maybe you should delete your PETSC_ARCH directory and recompile it?  I
> tested my branch. It should not that easily fail :)
>
> --Junchao Zhang
>
>
> On Tue, Mar 12, 2019 at 8:20 PM Manuel Valera  wrote:
>
>> Hi Mr Zhang, thanks for your reply,
>>
>> I just checked your branch out, reconfigured and recompiled and i am
>> still getting the same error from my last email (null argument, when
>> expected a valid pointer), do you have any idea why this can be happening?
>>
>> Thanks so much,
>>
>> Manuel
>>
>> On Tue, Mar 12, 2019 at 6:09 PM Zhang, Junchao 
>> wrote:
>>
>>> Manuel,
>>> I was working on a branch to revert the VecScatterCreate to
>>> VecScatterCreateWithData change. The change broke PETSc API and I think we
>>> do not need it. I had planed to do a pull request after my another PR is
>>> merged.
>>> But since it already affects you,  you can try this branch now, which is
>>> jczhang/fix-vecscattercreate-api
>>>
>>> Thanks.
>>> --Junchao Zhang
>>>
>>>
>>> On Tue, Mar 12, 2019 at 5:58 PM Jed Brown via petsc-users <
>>> petsc-users@mcs.anl.gov> wrote:
>>>
 Did you just update to 'master'?  See VecScatter changes:

 https://www.mcs.anl.gov/petsc/documentation/changes/dev.html

 Manuel Valera via petsc-users  writes:

 > Hello,
 >
 > I just updated petsc from the repo to the latest master branch
 version, and
 > a compilation problem popped up, it seems like the variable types are
 not
 > being acknowledged properly, what i have in a minimum working example
 > fashion is:
 >
 > #include 
 >> #include 
 >> #include 
 >> #include 
 >> #include 
 >> USE petscvec
 >> USE petscdmda
 >> USE petscdm
 >> USE petscis
 >> USE petscksp
 >> IS :: ScalarIS
 >> IS :: DummyIS
 >> VecScatter :: LargerToSmaller,to0,from0
 >> VecScatter :: SmallerToLarger
 >> PetscInt, ALLOCATABLE  :: pScalarDA(:), pDummyDA(:)
 >> PetscScalar:: rtol
 >> Vec:: Vec1
 >> Vec:: Vec2
 >> ! Create index sets
 >> allocate( pScalarDA(0:(gridx-1)*(gridy-1)*(gridz-1)-1) ,
 >> pDummyDA(0:(gridx-1)*(gridy-1)*(gridz-1)-1) )
 >> iter=0
 >> do k=0,gridz-2
 >> kplane = k*gridx*gridy
 >> do j=0,gridy-2
 >> do i=0,gridx-2
 >> pScalarDA(iter) = kplane + j*(gridx) + i
 >> iter = iter+1
 >> enddo
 >> enddo
 >> enddo
 >> pDummyDA = (/ (ind,
 ind=0,((gridx-1)*(gridy-1)*(gridz-1))-1) /)
 >> call
 >> ISCreateGeneral(PETSC_COMM_WORLD,(gridx-1)*(gridy-1)*(gridz-1), &
 >>
 >>  pScalarDA,PETSC_COPY_VALUES,ScalarIS,ierr)
 >> call
 >> ISCreateGeneral(PETSC_COMM_WORLD,(gridx-1)*(gridy-1)*(gridz-1), &
 >>
 >>  pDummyDA,PETSC_COPY_VALUES,DummyIS,ierr)
 >> deallocate(pScalarDA,pDummyDA, STAT=ierr)
 >> ! Create VecScatter contexts: LargerToSmaller &
 SmallerToLarger
 >> call DMDACreateNaturalVector(daScalars,Vec1,ierr)
 >> call DMDACreateNaturalVector(daDummy,Vec2,ierr)
 >> call
 >> VecScatterCreate(Vec1,ScalarIS,Vec2,DummyIS,LargerToSmaller,ierr)
 >> call
 >> VecScatterCreate(Vec2,DummyIS,Vec1,ScalarIS,SmallerToLarger,ierr)
 >> call VecDestroy(Vec1,ierr)
 >> call VecDestroy(Vec2,ierr)
 >
 >
 > And the error i get is the part i cannot really understand:
 >
 > matrixobjs.f90:99.34:
 >> call
 >> VecScatterCreate(Vec1,ScalarIS,Vec2,DummyIS,LargerToSmaller,ie
 >>  1
 >> Error: Type mismatch in argument 'a' at (1); passed TYPE(tvec) to
 >> INTEGER(4)
 >> matrixobjs.f90:100.34:
 >> call
 >> VecScatterCreate(Vec2,DummyIS,Vec1,ScalarIS,SmallerToLarger,ie
 >>  1
 >> Error: Type mismatch in argument 'a' at (1); passed TYPE(tvec) to
 >> INTEGER(4)
 >> make[1]: *** [matrixobjs.o] Error 1
 >> make[1]: Leaving directory `/usr/scratch/valera/ParGCCOM-Master/Src'
 >> make: *** [gcmSeamount] Error 2
 >
 >
 > What i find hard to understand is why/where my code is finding an
 integer
 > type? as you can see from the MWE header the variables types look
 correct,
 >
 > Any help is appreaciated,
 >
 > Thanks,

>>>


Re: [petsc-users] PetscScatterCreate type mismatch after update.

2019-03-12 Thread Zhang, Junchao via petsc-users
Maybe you should delete your PETSC_ARCH directory and recompile it?  I tested 
my branch. It should not that easily fail :)

--Junchao Zhang


On Tue, Mar 12, 2019 at 8:20 PM Manuel Valera 
mailto:mvaler...@sdsu.edu>> wrote:
Hi Mr Zhang, thanks for your reply,

I just checked your branch out, reconfigured and recompiled and i am still 
getting the same error from my last email (null argument, when expected a valid 
pointer), do you have any idea why this can be happening?

Thanks so much,

Manuel

On Tue, Mar 12, 2019 at 6:09 PM Zhang, Junchao 
mailto:jczh...@mcs.anl.gov>> wrote:
Manuel,
I was working on a branch to revert the VecScatterCreate to 
VecScatterCreateWithData change. The change broke PETSc API and I think we do 
not need it. I had planed to do a pull request after my another PR is merged.
But since it already affects you,  you can try this branch now, which is 
jczhang/fix-vecscattercreate-api

Thanks.
--Junchao Zhang


On Tue, Mar 12, 2019 at 5:58 PM Jed Brown via petsc-users 
mailto:petsc-users@mcs.anl.gov>> wrote:
Did you just update to 'master'?  See VecScatter changes:

https://www.mcs.anl.gov/petsc/documentation/changes/dev.html

Manuel Valera via petsc-users 
mailto:petsc-users@mcs.anl.gov>> writes:

> Hello,
>
> I just updated petsc from the repo to the latest master branch version, and
> a compilation problem popped up, it seems like the variable types are not
> being acknowledged properly, what i have in a minimum working example
> fashion is:
>
> #include 
>> #include 
>> #include 
>> #include 
>> #include 
>> USE petscvec
>> USE petscdmda
>> USE petscdm
>> USE petscis
>> USE petscksp
>> IS :: ScalarIS
>> IS :: DummyIS
>> VecScatter :: LargerToSmaller,to0,from0
>> VecScatter :: SmallerToLarger
>> PetscInt, ALLOCATABLE  :: pScalarDA(:), pDummyDA(:)
>> PetscScalar:: rtol
>> Vec:: Vec1
>> Vec:: Vec2
>> ! Create index sets
>> allocate( pScalarDA(0:(gridx-1)*(gridy-1)*(gridz-1)-1) ,
>> pDummyDA(0:(gridx-1)*(gridy-1)*(gridz-1)-1) )
>> iter=0
>> do k=0,gridz-2
>> kplane = k*gridx*gridy
>> do j=0,gridy-2
>> do i=0,gridx-2
>> pScalarDA(iter) = kplane + j*(gridx) + i
>> iter = iter+1
>> enddo
>> enddo
>> enddo
>> pDummyDA = (/ (ind, ind=0,((gridx-1)*(gridy-1)*(gridz-1))-1) /)
>> call
>> ISCreateGeneral(PETSC_COMM_WORLD,(gridx-1)*(gridy-1)*(gridz-1), &
>>
>>  pScalarDA,PETSC_COPY_VALUES,ScalarIS,ierr)
>> call
>> ISCreateGeneral(PETSC_COMM_WORLD,(gridx-1)*(gridy-1)*(gridz-1), &
>>
>>  pDummyDA,PETSC_COPY_VALUES,DummyIS,ierr)
>> deallocate(pScalarDA,pDummyDA, STAT=ierr)
>> ! Create VecScatter contexts: LargerToSmaller & SmallerToLarger
>> call DMDACreateNaturalVector(daScalars,Vec1,ierr)
>> call DMDACreateNaturalVector(daDummy,Vec2,ierr)
>> call
>> VecScatterCreate(Vec1,ScalarIS,Vec2,DummyIS,LargerToSmaller,ierr)
>> call
>> VecScatterCreate(Vec2,DummyIS,Vec1,ScalarIS,SmallerToLarger,ierr)
>> call VecDestroy(Vec1,ierr)
>> call VecDestroy(Vec2,ierr)
>
>
> And the error i get is the part i cannot really understand:
>
> matrixobjs.f90:99.34:
>> call
>> VecScatterCreate(Vec1,ScalarIS,Vec2,DummyIS,LargerToSmaller,ie
>>  1
>> Error: Type mismatch in argument 'a' at (1); passed TYPE(tvec) to
>> INTEGER(4)
>> matrixobjs.f90:100.34:
>> call
>> VecScatterCreate(Vec2,DummyIS,Vec1,ScalarIS,SmallerToLarger,ie
>>  1
>> Error: Type mismatch in argument 'a' at (1); passed TYPE(tvec) to
>> INTEGER(4)
>> make[1]: *** [matrixobjs.o] Error 1
>> make[1]: Leaving directory `/usr/scratch/valera/ParGCCOM-Master/Src'
>> make: *** [gcmSeamount] Error 2
>
>
> What i find hard to understand is why/where my code is finding an integer
> type? as you can see from the MWE header the variables types look correct,
>
> Any help is appreaciated,
>
> Thanks,


Re: [petsc-users] PetscScatterCreate type mismatch after update.

2019-03-12 Thread Manuel Valera via petsc-users
Hi Mr Zhang, thanks for your reply,

I just checked your branch out, reconfigured and recompiled and i am still
getting the same error from my last email (null argument, when expected a
valid pointer), do you have any idea why this can be happening?

Thanks so much,

Manuel

On Tue, Mar 12, 2019 at 6:09 PM Zhang, Junchao  wrote:

> Manuel,
> I was working on a branch to revert the VecScatterCreate to
> VecScatterCreateWithData change. The change broke PETSc API and I think we
> do not need it. I had planed to do a pull request after my another PR is
> merged.
> But since it already affects you,  you can try this branch now, which is
> jczhang/fix-vecscattercreate-api
>
> Thanks.
> --Junchao Zhang
>
>
> On Tue, Mar 12, 2019 at 5:58 PM Jed Brown via petsc-users <
> petsc-users@mcs.anl.gov> wrote:
>
>> Did you just update to 'master'?  See VecScatter changes:
>>
>> https://www.mcs.anl.gov/petsc/documentation/changes/dev.html
>>
>> Manuel Valera via petsc-users  writes:
>>
>> > Hello,
>> >
>> > I just updated petsc from the repo to the latest master branch version,
>> and
>> > a compilation problem popped up, it seems like the variable types are
>> not
>> > being acknowledged properly, what i have in a minimum working example
>> > fashion is:
>> >
>> > #include 
>> >> #include 
>> >> #include 
>> >> #include 
>> >> #include 
>> >> USE petscvec
>> >> USE petscdmda
>> >> USE petscdm
>> >> USE petscis
>> >> USE petscksp
>> >> IS :: ScalarIS
>> >> IS :: DummyIS
>> >> VecScatter :: LargerToSmaller,to0,from0
>> >> VecScatter :: SmallerToLarger
>> >> PetscInt, ALLOCATABLE  :: pScalarDA(:), pDummyDA(:)
>> >> PetscScalar:: rtol
>> >> Vec:: Vec1
>> >> Vec:: Vec2
>> >> ! Create index sets
>> >> allocate( pScalarDA(0:(gridx-1)*(gridy-1)*(gridz-1)-1) ,
>> >> pDummyDA(0:(gridx-1)*(gridy-1)*(gridz-1)-1) )
>> >> iter=0
>> >> do k=0,gridz-2
>> >> kplane = k*gridx*gridy
>> >> do j=0,gridy-2
>> >> do i=0,gridx-2
>> >> pScalarDA(iter) = kplane + j*(gridx) + i
>> >> iter = iter+1
>> >> enddo
>> >> enddo
>> >> enddo
>> >> pDummyDA = (/ (ind,
>> ind=0,((gridx-1)*(gridy-1)*(gridz-1))-1) /)
>> >> call
>> >> ISCreateGeneral(PETSC_COMM_WORLD,(gridx-1)*(gridy-1)*(gridz-1), &
>> >>
>> >>  pScalarDA,PETSC_COPY_VALUES,ScalarIS,ierr)
>> >> call
>> >> ISCreateGeneral(PETSC_COMM_WORLD,(gridx-1)*(gridy-1)*(gridz-1), &
>> >>
>> >>  pDummyDA,PETSC_COPY_VALUES,DummyIS,ierr)
>> >> deallocate(pScalarDA,pDummyDA, STAT=ierr)
>> >> ! Create VecScatter contexts: LargerToSmaller &
>> SmallerToLarger
>> >> call DMDACreateNaturalVector(daScalars,Vec1,ierr)
>> >> call DMDACreateNaturalVector(daDummy,Vec2,ierr)
>> >> call
>> >> VecScatterCreate(Vec1,ScalarIS,Vec2,DummyIS,LargerToSmaller,ierr)
>> >> call
>> >> VecScatterCreate(Vec2,DummyIS,Vec1,ScalarIS,SmallerToLarger,ierr)
>> >> call VecDestroy(Vec1,ierr)
>> >> call VecDestroy(Vec2,ierr)
>> >
>> >
>> > And the error i get is the part i cannot really understand:
>> >
>> > matrixobjs.f90:99.34:
>> >> call
>> >> VecScatterCreate(Vec1,ScalarIS,Vec2,DummyIS,LargerToSmaller,ie
>> >>  1
>> >> Error: Type mismatch in argument 'a' at (1); passed TYPE(tvec) to
>> >> INTEGER(4)
>> >> matrixobjs.f90:100.34:
>> >> call
>> >> VecScatterCreate(Vec2,DummyIS,Vec1,ScalarIS,SmallerToLarger,ie
>> >>  1
>> >> Error: Type mismatch in argument 'a' at (1); passed TYPE(tvec) to
>> >> INTEGER(4)
>> >> make[1]: *** [matrixobjs.o] Error 1
>> >> make[1]: Leaving directory `/usr/scratch/valera/ParGCCOM-Master/Src'
>> >> make: *** [gcmSeamount] Error 2
>> >
>> >
>> > What i find hard to understand is why/where my code is finding an
>> integer
>> > type? as you can see from the MWE header the variables types look
>> correct,
>> >
>> > Any help is appreaciated,
>> >
>> > Thanks,
>>
>


Re: [petsc-users] PetscScatterCreate type mismatch after update.

2019-03-12 Thread Zhang, Junchao via petsc-users
Manuel,
I was working on a branch to revert the VecScatterCreate to 
VecScatterCreateWithData change. The change broke PETSc API and I think we do 
not need it. I had planed to do a pull request after my another PR is merged.
But since it already affects you,  you can try this branch now, which is 
jczhang/fix-vecscattercreate-api

Thanks.
--Junchao Zhang


On Tue, Mar 12, 2019 at 5:58 PM Jed Brown via petsc-users 
mailto:petsc-users@mcs.anl.gov>> wrote:
Did you just update to 'master'?  See VecScatter changes:

https://www.mcs.anl.gov/petsc/documentation/changes/dev.html

Manuel Valera via petsc-users 
mailto:petsc-users@mcs.anl.gov>> writes:

> Hello,
>
> I just updated petsc from the repo to the latest master branch version, and
> a compilation problem popped up, it seems like the variable types are not
> being acknowledged properly, what i have in a minimum working example
> fashion is:
>
> #include 
>> #include 
>> #include 
>> #include 
>> #include 
>> USE petscvec
>> USE petscdmda
>> USE petscdm
>> USE petscis
>> USE petscksp
>> IS :: ScalarIS
>> IS :: DummyIS
>> VecScatter :: LargerToSmaller,to0,from0
>> VecScatter :: SmallerToLarger
>> PetscInt, ALLOCATABLE  :: pScalarDA(:), pDummyDA(:)
>> PetscScalar:: rtol
>> Vec:: Vec1
>> Vec:: Vec2
>> ! Create index sets
>> allocate( pScalarDA(0:(gridx-1)*(gridy-1)*(gridz-1)-1) ,
>> pDummyDA(0:(gridx-1)*(gridy-1)*(gridz-1)-1) )
>> iter=0
>> do k=0,gridz-2
>> kplane = k*gridx*gridy
>> do j=0,gridy-2
>> do i=0,gridx-2
>> pScalarDA(iter) = kplane + j*(gridx) + i
>> iter = iter+1
>> enddo
>> enddo
>> enddo
>> pDummyDA = (/ (ind, ind=0,((gridx-1)*(gridy-1)*(gridz-1))-1) /)
>> call
>> ISCreateGeneral(PETSC_COMM_WORLD,(gridx-1)*(gridy-1)*(gridz-1), &
>>
>>  pScalarDA,PETSC_COPY_VALUES,ScalarIS,ierr)
>> call
>> ISCreateGeneral(PETSC_COMM_WORLD,(gridx-1)*(gridy-1)*(gridz-1), &
>>
>>  pDummyDA,PETSC_COPY_VALUES,DummyIS,ierr)
>> deallocate(pScalarDA,pDummyDA, STAT=ierr)
>> ! Create VecScatter contexts: LargerToSmaller & SmallerToLarger
>> call DMDACreateNaturalVector(daScalars,Vec1,ierr)
>> call DMDACreateNaturalVector(daDummy,Vec2,ierr)
>> call
>> VecScatterCreate(Vec1,ScalarIS,Vec2,DummyIS,LargerToSmaller,ierr)
>> call
>> VecScatterCreate(Vec2,DummyIS,Vec1,ScalarIS,SmallerToLarger,ierr)
>> call VecDestroy(Vec1,ierr)
>> call VecDestroy(Vec2,ierr)
>
>
> And the error i get is the part i cannot really understand:
>
> matrixobjs.f90:99.34:
>> call
>> VecScatterCreate(Vec1,ScalarIS,Vec2,DummyIS,LargerToSmaller,ie
>>  1
>> Error: Type mismatch in argument 'a' at (1); passed TYPE(tvec) to
>> INTEGER(4)
>> matrixobjs.f90:100.34:
>> call
>> VecScatterCreate(Vec2,DummyIS,Vec1,ScalarIS,SmallerToLarger,ie
>>  1
>> Error: Type mismatch in argument 'a' at (1); passed TYPE(tvec) to
>> INTEGER(4)
>> make[1]: *** [matrixobjs.o] Error 1
>> make[1]: Leaving directory `/usr/scratch/valera/ParGCCOM-Master/Src'
>> make: *** [gcmSeamount] Error 2
>
>
> What i find hard to understand is why/where my code is finding an integer
> type? as you can see from the MWE header the variables types look correct,
>
> Any help is appreaciated,
>
> Thanks,


Re: [petsc-users] PetscScatterCreate type mismatch after update.

2019-03-12 Thread Manuel Valera via petsc-users
Hi,

So i just solved that problem but now it looks my code broke somewhere
else, i have a script in place to scatter/gather the information to root in
order to write it to a file (i know, we need to make this parallel I/O but
that's future work). Such script looks like this:


>   SUBROUTINE WriteToFile_grid()
>
> PetscErrorCode:: ierrp
> PetscMPIInt   :: rank
> PetscInt  :: iter
> Vec   ::
> CenterX,CenterY,CenterZ,Nat1,Nat2,seqvec
> PetscScalar, pointer  :: tmp3d(:,:,:),tmp4d(:,:,:,:),arr(:)
> VecScatter:: LargerToSmaller, scatterctx
> INTEGER   :: i,j,k, ierr
>
> call MPI_Comm_rank(PETSC_COMM_WORLD, rank, ierrp)
>
> !#!
> ! Grid Cell Centers: x-component
> !
>
> !#!
> ! Extract x-component
> call DMDAVecGetArrayF90(daGrid,GridCenters,tmp4d,ierrp)
> call DMCreateGlobalVector(daSingle,CenterX,ierrp)
> call DMDAVecGetArrayF90(daSingle,CenterX,tmp3d,ierrp)
> tmp3d(:,:,:) = tmp4d(0,:,:,:)
> call DMDAVecRestoreArrayF90(daSingle,CenterX,tmp3d,ierrp)
> call DMDAVecRestoreArrayF90(daGrid,GridCenters,tmp4d,ierrp)
> ! Scatter to daWriteCenters
> call DMDACreateNaturalVector(daSingle,Nat1,ierrp)
> call
> DMDAGlobalToNaturalBegin(daSingle,CenterX,INSERT_VALUES,Nat1,ierrp)
> call
> DMDAGlobalToNaturalEnd(daSingle,CenterX,INSERT_VALUES,Nat1,ierrp)
> call VecDestroy(CenterX,ierrp)
> call DMDACreateNaturalVector(daWriteCenters,Nat2,ierrp)
> call
> VecScatterCreate(Nat1,SingleIS,Nat2,WriteIS,LargerToSmaller,ierrp)
> call
> VecScatterBegin(LargerToSmaller,Nat1,Nat2,INSERT_VALUES,SCATTER_FORWARD,ierrp)
> call
> VecScatterEnd(LargerToSmaller,Nat1,Nat2,INSERT_VALUES,SCATTER_FORWARD,ierrp)
> call VecScatterDestroy(LargerToSmaller,ierrp)
> call VecDestroy(Nat1,ierrp)
> ! Send to root
> call VecScatterCreateToZero(Nat2,scatterctx,seqvec,ierrp)
> call
> VecScatterBegin(scatterctx,Nat2,seqvec,INSERT_VALUES,SCATTER_FORWARD,ierrp)
> call
> VecScatterEnd(scatterctx,Nat2,seqvec,INSERT_VALUES,SCATTER_FORWARD,ierrp)
> call VecScatterDestroy(scatterctx,ierrp)
>call VecDestroy(Nat2,ierrp)
> ! Let root write to netCDF file
> if (rank == 0) then
> allocate(buffer(1:IMax-1,1:JMax-1,1:KMax-1),STAT=ierr)
> call VecGetArrayReadF90(seqvec,arr,ierrp)
> iter = 1
> do k=1,KMax-1
> do j=1,JMax-1
> do i=1,IMax-1
> buffer(i,j,k) = arr(iter)
> iter = iter + 1
> enddo
> enddo
> enddo
> call VecRestoreArrayReadF90(seqvec,arr,ierrp)
> call
> nc_check(nf90_put_var(ncid,xID,buffer,start=(/1,1,1/),count=(/IMax-1,JMax-1,KMax-1/)),
> &
>   'WriteNetCDF', context='put_var GridCenterX
> in '//trim(output_filename))
> deallocate(buffer,STAT=ierr)
> endif
> call VecDestroy(seqvec,ierrp)


And then the process is repeated for each variable to output. Notice the
vector seqvec is being destroyed at the end. Using petsc v3.10.0 this
script worked without problems. After updating to v3.10.4 it no longer
works. Gives the following error:

[0]PETSC ERROR: - Error Message
> --
> [0]PETSC ERROR: Null argument, when expecting valid pointer
> [0]PETSC ERROR: Null Object: Parameter # 3
> [0]PETSC ERROR: See http://www.mcs.anl.gov/petsc/documentation/faq.html
> for trouble shooting.
> [0]PETSC ERROR: Petsc Release Version 3.10.4, unknown
> [0]PETSC ERROR: ./gcmSeamount on a petsc-debug named ocean by valera Tue
> Mar 12 17:59:43 2019
> [0]PETSC ERROR: Configure options --known-level1-dcache-size=32768
> --known-level1-dcache-linesize=64 --known-level1-dcache-assoc=8
> --known-sizeof-char=1 --known-sizeof-void-p=8 --known-sizeof-short=2
> --known-sizeof-int=4 --known-sizeof-long=8 --known-sizeof-long-long=8
> --known-sizeof-float=4 --known-sizeof-double=8 --known-sizeof-size_t=8
> --known-bits-per-byte=8 --known-memcmp-ok=1 --known-sizeof-MPI_Comm=8
> --known-sizeof-MPI_Fint=4 --known-mpi-long-double=1 --known-mpi-int64_t=1
> --known-mpi-c-double-complex=1 --known-has-attribute-aligned=1
> PETSC_ARCH=petsc-debug --COPTFLAGS=-O2 --CXXOPTFLAGS=-O2 --FOPTFLAGS=-O2
> --with-cc=mpicc 

Re: [petsc-users] PetscScatterCreate type mismatch after update.

2019-03-12 Thread Jed Brown via petsc-users
Did you just update to 'master'?  See VecScatter changes:

https://www.mcs.anl.gov/petsc/documentation/changes/dev.html

Manuel Valera via petsc-users  writes:

> Hello,
>
> I just updated petsc from the repo to the latest master branch version, and
> a compilation problem popped up, it seems like the variable types are not
> being acknowledged properly, what i have in a minimum working example
> fashion is:
>
> #include 
>> #include 
>> #include 
>> #include 
>> #include 
>> USE petscvec
>> USE petscdmda
>> USE petscdm
>> USE petscis
>> USE petscksp
>> IS :: ScalarIS
>> IS :: DummyIS
>> VecScatter :: LargerToSmaller,to0,from0
>> VecScatter :: SmallerToLarger
>> PetscInt, ALLOCATABLE  :: pScalarDA(:), pDummyDA(:)
>> PetscScalar:: rtol
>> Vec:: Vec1
>> Vec:: Vec2
>> ! Create index sets
>> allocate( pScalarDA(0:(gridx-1)*(gridy-1)*(gridz-1)-1) ,
>> pDummyDA(0:(gridx-1)*(gridy-1)*(gridz-1)-1) )
>> iter=0
>> do k=0,gridz-2
>> kplane = k*gridx*gridy
>> do j=0,gridy-2
>> do i=0,gridx-2
>> pScalarDA(iter) = kplane + j*(gridx) + i
>> iter = iter+1
>> enddo
>> enddo
>> enddo
>> pDummyDA = (/ (ind, ind=0,((gridx-1)*(gridy-1)*(gridz-1))-1) /)
>> call
>> ISCreateGeneral(PETSC_COMM_WORLD,(gridx-1)*(gridy-1)*(gridz-1), &
>>
>>  pScalarDA,PETSC_COPY_VALUES,ScalarIS,ierr)
>> call
>> ISCreateGeneral(PETSC_COMM_WORLD,(gridx-1)*(gridy-1)*(gridz-1), &
>>
>>  pDummyDA,PETSC_COPY_VALUES,DummyIS,ierr)
>> deallocate(pScalarDA,pDummyDA, STAT=ierr)
>> ! Create VecScatter contexts: LargerToSmaller & SmallerToLarger
>> call DMDACreateNaturalVector(daScalars,Vec1,ierr)
>> call DMDACreateNaturalVector(daDummy,Vec2,ierr)
>> call
>> VecScatterCreate(Vec1,ScalarIS,Vec2,DummyIS,LargerToSmaller,ierr)
>> call
>> VecScatterCreate(Vec2,DummyIS,Vec1,ScalarIS,SmallerToLarger,ierr)
>> call VecDestroy(Vec1,ierr)
>> call VecDestroy(Vec2,ierr)
>
>
> And the error i get is the part i cannot really understand:
>
> matrixobjs.f90:99.34:
>> call
>> VecScatterCreate(Vec1,ScalarIS,Vec2,DummyIS,LargerToSmaller,ie
>>  1
>> Error: Type mismatch in argument 'a' at (1); passed TYPE(tvec) to
>> INTEGER(4)
>> matrixobjs.f90:100.34:
>> call
>> VecScatterCreate(Vec2,DummyIS,Vec1,ScalarIS,SmallerToLarger,ie
>>  1
>> Error: Type mismatch in argument 'a' at (1); passed TYPE(tvec) to
>> INTEGER(4)
>> make[1]: *** [matrixobjs.o] Error 1
>> make[1]: Leaving directory `/usr/scratch/valera/ParGCCOM-Master/Src'
>> make: *** [gcmSeamount] Error 2
>
>
> What i find hard to understand is why/where my code is finding an integer
> type? as you can see from the MWE header the variables types look correct,
>
> Any help is appreaciated,
>
> Thanks,


[petsc-users] PetscScatterCreate type mismatch after update.

2019-03-12 Thread Manuel Valera via petsc-users
Hello,

I just updated petsc from the repo to the latest master branch version, and
a compilation problem popped up, it seems like the variable types are not
being acknowledged properly, what i have in a minimum working example
fashion is:

#include 
> #include 
> #include 
> #include 
> #include 
> USE petscvec
> USE petscdmda
> USE petscdm
> USE petscis
> USE petscksp
> IS :: ScalarIS
> IS :: DummyIS
> VecScatter :: LargerToSmaller,to0,from0
> VecScatter :: SmallerToLarger
> PetscInt, ALLOCATABLE  :: pScalarDA(:), pDummyDA(:)
> PetscScalar:: rtol
> Vec:: Vec1
> Vec:: Vec2
> ! Create index sets
> allocate( pScalarDA(0:(gridx-1)*(gridy-1)*(gridz-1)-1) ,
> pDummyDA(0:(gridx-1)*(gridy-1)*(gridz-1)-1) )
> iter=0
> do k=0,gridz-2
> kplane = k*gridx*gridy
> do j=0,gridy-2
> do i=0,gridx-2
> pScalarDA(iter) = kplane + j*(gridx) + i
> iter = iter+1
> enddo
> enddo
> enddo
> pDummyDA = (/ (ind, ind=0,((gridx-1)*(gridy-1)*(gridz-1))-1) /)
> call
> ISCreateGeneral(PETSC_COMM_WORLD,(gridx-1)*(gridy-1)*(gridz-1), &
>
>  pScalarDA,PETSC_COPY_VALUES,ScalarIS,ierr)
> call
> ISCreateGeneral(PETSC_COMM_WORLD,(gridx-1)*(gridy-1)*(gridz-1), &
>
>  pDummyDA,PETSC_COPY_VALUES,DummyIS,ierr)
> deallocate(pScalarDA,pDummyDA, STAT=ierr)
> ! Create VecScatter contexts: LargerToSmaller & SmallerToLarger
> call DMDACreateNaturalVector(daScalars,Vec1,ierr)
> call DMDACreateNaturalVector(daDummy,Vec2,ierr)
> call
> VecScatterCreate(Vec1,ScalarIS,Vec2,DummyIS,LargerToSmaller,ierr)
> call
> VecScatterCreate(Vec2,DummyIS,Vec1,ScalarIS,SmallerToLarger,ierr)
> call VecDestroy(Vec1,ierr)
> call VecDestroy(Vec2,ierr)


And the error i get is the part i cannot really understand:

matrixobjs.f90:99.34:
> call
> VecScatterCreate(Vec1,ScalarIS,Vec2,DummyIS,LargerToSmaller,ie
>  1
> Error: Type mismatch in argument 'a' at (1); passed TYPE(tvec) to
> INTEGER(4)
> matrixobjs.f90:100.34:
> call
> VecScatterCreate(Vec2,DummyIS,Vec1,ScalarIS,SmallerToLarger,ie
>  1
> Error: Type mismatch in argument 'a' at (1); passed TYPE(tvec) to
> INTEGER(4)
> make[1]: *** [matrixobjs.o] Error 1
> make[1]: Leaving directory `/usr/scratch/valera/ParGCCOM-Master/Src'
> make: *** [gcmSeamount] Error 2


What i find hard to understand is why/where my code is finding an integer
type? as you can see from the MWE header the variables types look correct,

Any help is appreaciated,

Thanks,


Re: [petsc-users] PCFieldSplit with MatNest

2019-03-12 Thread Zhang, Junchao via petsc-users
Hi, Manuel,
  I recently fixed a problem in VecRestoreArrayRead. Basically, I added 
VecRestoreArrayRead_Nest. Could you try the master branch of PETSc to see if it 
fixes your problem?
  Thanks.

--Junchao Zhang


On Mon, Mar 11, 2019 at 6:56 AM Manuel Colera Rico via petsc-users 
mailto:petsc-users@mcs.anl.gov>> wrote:
Hello,

I need to solve a 2*2 block linear system. The matrices A_00, A_01,
A_10, A_11 are constructed separately via MatCreateSeqAIJWithArrays and
MatCreateSeqSBAIJWithArrays. Then, I construct the full system matrix
with MatCreateNest, and use MatNestGetISs and PCFieldSplitSetIS to set
up the PC, trying to follow the procedure described here:
https://www.mcs.anl.gov/petsc/petsc-current/src/snes/examples/tutorials/ex70.c.html.

However, when I run the code with Leak Sanitizer, I get the following error:

=
==54927==ERROR: AddressSanitizer: attempting free on address which was
not malloc()-ed: 0x62751ab8 in thread T0
 #0 0x7fbd95c08f30 in __interceptor_free
../../../../gcc-8.1.0/libsanitizer/asan/asan_malloc_linux.cc:66
 #1 0x7fbd92b99dcd in PetscFreeAlign
(/opt/PETSc_library/petsc/manuel_OpenBLAS_petsc/lib/libpetsc.so.3.8+0x146dcd)
 #2 0x7fbd92ce0178 in VecRestoreArray_Nest
(/opt/PETSc_library/petsc/manuel_OpenBLAS_petsc/lib/libpetsc.so.3.8+0x28d178)
 #3 0x7fbd92cd627d in VecRestoreArrayRead
(/opt/PETSc_library/petsc/manuel_OpenBLAS_petsc/lib/libpetsc.so.3.8+0x28327d)
 #4 0x7fbd92d1189e in VecScatterBegin_SSToSS
(/opt/PETSc_library/petsc/manuel_OpenBLAS_petsc/lib/libpetsc.so.3.8+0x2be89e)
 #5 0x7fbd92d1a414 in VecScatterBegin
(/opt/PETSc_library/petsc/manuel_OpenBLAS_petsc/lib/libpetsc.so.3.8+0x2c7414)
 #6 0x7fbd934a999c in PCApply_FieldSplit
(/opt/PETSc_library/petsc/manuel_OpenBLAS_petsc/lib/libpetsc.so.3.8+0xa5699c)
 #7 0x7fbd93369071 in PCApply
(/opt/PETSc_library/petsc/manuel_OpenBLAS_petsc/lib/libpetsc.so.3.8+0x916071)
 #8 0x7fbd934efe77 in KSPInitialResidual
(/opt/PETSc_library/petsc/manuel_OpenBLAS_petsc/lib/libpetsc.so.3.8+0xa9ce77)
 #9 0x7fbd9350272c in KSPSolve_GMRES
(/opt/PETSc_library/petsc/manuel_OpenBLAS_petsc/lib/libpetsc.so.3.8+0xaaf72c)
 #10 0x7fbd934e3c01 in KSPSolve
(/opt/PETSc_library/petsc/manuel_OpenBLAS_petsc/lib/libpetsc.so.3.8+0xa90c01)

Disabling Leak Sanitizer also outputs an "invalid pointer" error.

Did I forget something when writing the code?

Thank you,

Manuel

---



Re: [petsc-users] Matrix Partitioning using PARMETIS

2019-03-12 Thread Smith, Barry F. via petsc-users


   Yes, by default there is one subdomain per process so if you run on one 
process you will get all zero indices. Run on two processes and you should see 
a partitioning.  See also MatPartitioningSetNParts()

   Barry


> On Mar 12, 2019, at 3:03 AM, Eda Oktay via petsc-users 
>  wrote:
> 
> Hello,
> 
> I have a Laplacian matrix PL of matrix A and I try to partition A using 
> PARMETIS. Since PL is sequential and not adjacency matrix, I converted PL to 
> AL, then write the following code:
> 
>   ierr = MatConvert(PL,MATMPIADJ,MAT_INITIAL_MATRIX,);CHKERRQ(ierr);   
>  
>   ierr = MatMeshToCellGraph(AL,2,);CHKERRQ(ierr); 
>   ierr = MatPartitioningCreate(MPI_COMM_WORLD,);CHKERRQ(ierr);
>   ierr = MatPartitioningSetAdjacency(part,dual);CHKERRQ(ierr);   
>   ierr = MatPartitioningSetFromOptions(part);CHKERRQ(ierr);
>   ierr = MatPartitioningApply(part,);CHKERRQ(ierr);
>   ierr = ISView(partitioning,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);   
>   ierr = ISDestroy();CHKERRQ(ierr);
>   ierr = MatPartitioningDestroy();CHKERRQ(ierr);
> 
> However, when I look at partitioning with ISView, the index set consists of 
> zeros only. Is that because I have only one processor and my codes are 
> written for only one processor, or is there another problem? I ran my code 
> with -mat_partitioning_type parmetis.
> 
> Thanks,
> 
> Eda



Re: [petsc-users] DMPlexComputeCellGeometryFVM: "Cannot handle faces with 1 vertices"

2019-03-12 Thread Chris Finn via petsc-users
I think I found a hint in src/dm/impls/plex/examples/tests/ex5.c to get my code 
corrected.
Apparently the contiguous numbering is
NOT: cells, vertices, edges, faces (as you said),
BUT: cells, vertices, faces, edges.

Then in numPoints, you don't put depth 0..3 or depth 3..0,
BUT: #vertices, #faces, #edges, #cells.

This seems completely random, but works in the attached example.
Maybe you should finally take the time to document this mess if even the 
developers no longer understand it.

regards
Chris

Mar 12, 2019, 1:38 PM by knep...@gmail.com:

> On Tue, Mar 12, 2019 at 9:25 AM <> finnkochin...@keemail.me 
> > > wrote:
>
>> Thank you,
>> I added DMPlexInterpolate(). Now I get a different error in 
>> DMPlexInterpolate():
>> ...
>>
>
> Okay, you did not need interpolate. You already specified all the levels. 
> However, your orientations are wrong.
> Get rid of that code.
>  
>
>>
>> [0]PETSC ERROR: Dimension 0 not supported
>> ...
>>
>> The source code and full output are attached. Anybody able to fix this?
>>
>
> The problem is a numbering convention in the library. Plex will accept any 
> consistent DAG. However, if you
> want to use other things in the library, like geometry routines, then there 
> is a convention on numbering (which
> makes many things simpler). We require that you number contiguously:
>
>   Cells:      [0, Nc)
>   Vertices: [Nc, Nc+Nv)
>   Edges:    [Nc+Nv, Nc+Nv+Ne)
>   Faces:    [Nc+Nv+Ne, Nc+Nv+Ne+Nf)
>
> You numbered the faces in the vertices slot,  so the geometry routines got 
> confused.
>
>   Thanks,
>
>     Matt
>  
>
>>
>> regards
>> Chris
>>
>> -- 
>> Securely sent with Tutanota. Get your own encrypted, ad-free mailbox: 
>> https://tutanota.com 
>>
>>
>> Mar 12, 2019, 1:08 PM by >> knep...@gmail.com >> :
>>
>>> On Tue, Mar 12, 2019 at 9:03 AM Chris Finn via petsc-users <>>> 
>>> petsc-users@mcs.anl.gov  wrote:
>>>
 Hello,
 with the code below, I create a tetrahedron using DMPlexCreateFromDAG, 
 then I try to run DMPlexComputeCellGeometryFVM on this cell. The latter 
 call fails with output:

>>>
>>> All the geometry stuff requires that you interpolate the mesh I think, Just 
>>> use DMPlexInterpolate().
>>>
>>>   Thanks,
>>>
>>>     Matt
>>>  
>>>

 ...
 [0]PETSC ERROR: - Error Message -
 [0]PETSC ERROR: Argument out of range
 [0]PETSC ERROR: Cannot handle faces with 1 vertices
 ...
 (full output is attached).

 What am I doing wrong?
 regards
 Chris

 Here is the code:
 static char help[] = "No help \n";

 #include 

 #undef __FUNCT__
 #define __FUNCT__ "main"
 int main(int argc,char **args)
 {
   PetscErrorCode ierr;
   PetscInitialize(,,(char*)0,help);

   DM    dm;
   int cStart,cEnd;
   int fStart,fEnd;

   int depth = 3;
   int dim = 3;

   PetscInt    numPoints[4]    = {1,4,6,4};
   PetscInt    coneSize[15] = {4,3,3,3,3,2,2,2,2,2,2,0,0,0,0};
   PetscInt    cones[28]    = {1,2,3,4, 5,9,8, 9,6,10, 10,8,7, 
 5,6,7, 11,12, 12,13, 13,11, 11,14, 12,14, 13,14};
   PetscInt    coneOrientations[28] = {0 };
   PetscScalar vertexCoords[12] = {0,0,0, 1,0,0, 0,1,0, 0,0,1};

   DMCreate(PETSC_COMM_WORLD, );
   DMSetType(dm, DMPLEX);
   DMSetDimension(dm,dim);
   DMPlexCreateFromDAG(dm, depth, numPoints, coneSize, cones, 
 coneOrientations, vertexCoords);

   DMPlexGetHeightStratum(dm, 0, , );
   DMPlexGetHeightStratum(dm, 1, , );

   for (int k =cStart;k>>>     double vol;
     double centroid[3];
     double normal[3];
     ierr = DMPlexComputeCellGeometryFVM(dm, k, , 
 centroid,NULL);CHKERRQ(ierr);
     printf("FVM: V=%f c=(%f %f %f) n=(%f %f 
 %f)\n",vol,centroid[0],centroid[1],centroid[2],
   normal[0],normal[1],normal[2]);
   }
   ierr = PetscFinalize();
   return 0;
 }



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

static char help[] = "No help \n";

#include 

#undef __FUNCT__
#define __FUNCT__ "main"
int main(int argc,char **args)
{
  PetscErrorCode ierr;
  PetscInitialize(,,(char*)0,help);

  DMdm;
  int cStart,cEnd;
  int fStart,fEnd;

  int depth = 3;
  int dim = 3;


Re: [petsc-users] DMPlexComputeCellGeometryFVM: "Cannot handle faces with 1 vertices"

2019-03-12 Thread Chris Finn via petsc-users
Next try:
Vertex and face numbering corrected (I hope), code and full output attached.
Orientations still wrong, but the crash is elsewhere:
...
reconstructed depth,dim: 1 3
[0]PETSC ERROR: - Error Message 
--
[0]PETSC ERROR: Invalid argument
[0]PETSC ERROR: Mesh must be interpolated
...

DMPlexComputeCellGeometryFVM() checks depth and dim in the beginning, finds 
they are different and quits because mesh seems not interpolated.
Actually, DMPlexGetDepth and DMGetDimension return 1 and 3 with this 
vertex/face numbering. They returned the correct 3 and 3 with my previous wrong 
numbering.

It can't be so difficult to create a simplex from scratch?

regards
Chris

Mar 12, 2019, 1:38 PM by knep...@gmail.com:

> On Tue, Mar 12, 2019 at 9:25 AM <> finnkochin...@keemail.me 
> > > wrote:
>
>> Thank you,
>> I added DMPlexInterpolate(). Now I get a different error in 
>> DMPlexInterpolate():
>> ...
>>
>
> Okay, you did not need interpolate. You already specified all the levels. 
> However, your orientations are wrong.
> Get rid of that code.
>  
>
>>
>> [0]PETSC ERROR: Dimension 0 not supported
>> ...
>>
>> The source code and full output are attached. Anybody able to fix this?
>>
>
> The problem is a numbering convention in the library. Plex will accept any 
> consistent DAG. However, if you
> want to use other things in the library, like geometry routines, then there 
> is a convention on numbering (which
> makes many things simpler). We require that you number contiguously:
>
>   Cells:      [0, Nc)
>   Vertices: [Nc, Nc+Nv)
>   Edges:    [Nc+Nv, Nc+Nv+Ne)
>   Faces:    [Nc+Nv+Ne, Nc+Nv+Ne+Nf)
>
> You numbered the faces in the vertices slot,  so the geometry routines got 
> confused.
>
>   Thanks,
>
>     Matt
>  
>
>>
>> regards
>> Chris
>>
>> -- 
>> Securely sent with Tutanota. Get your own encrypted, ad-free mailbox: 
>> https://tutanota.com 
>>
>>
>> Mar 12, 2019, 1:08 PM by >> knep...@gmail.com >> :
>>
>>> On Tue, Mar 12, 2019 at 9:03 AM Chris Finn via petsc-users <>>> 
>>> petsc-users@mcs.anl.gov  wrote:
>>>
 Hello,
 with the code below, I create a tetrahedron using DMPlexCreateFromDAG, 
 then I try to run DMPlexComputeCellGeometryFVM on this cell. The latter 
 call fails with output:

>>>
>>> All the geometry stuff requires that you interpolate the mesh I think, Just 
>>> use DMPlexInterpolate().
>>>
>>>   Thanks,
>>>
>>>     Matt
>>>  
>>>

 ...
 [0]PETSC ERROR: - Error Message -
 [0]PETSC ERROR: Argument out of range
 [0]PETSC ERROR: Cannot handle faces with 1 vertices
 ...
 (full output is attached).

 What am I doing wrong?
 regards
 Chris

 Here is the code:
 static char help[] = "No help \n";

 #include 

 #undef __FUNCT__
 #define __FUNCT__ "main"
 int main(int argc,char **args)
 {
   PetscErrorCode ierr;
   PetscInitialize(,,(char*)0,help);

   DM    dm;
   int cStart,cEnd;
   int fStart,fEnd;

   int depth = 3;
   int dim = 3;

   PetscInt    numPoints[4]    = {1,4,6,4};
   PetscInt    coneSize[15] = {4,3,3,3,3,2,2,2,2,2,2,0,0,0,0};
   PetscInt    cones[28]    = {1,2,3,4, 5,9,8, 9,6,10, 10,8,7, 
 5,6,7, 11,12, 12,13, 13,11, 11,14, 12,14, 13,14};
   PetscInt    coneOrientations[28] = {0 };
   PetscScalar vertexCoords[12] = {0,0,0, 1,0,0, 0,1,0, 0,0,1};

   DMCreate(PETSC_COMM_WORLD, );
   DMSetType(dm, DMPLEX);
   DMSetDimension(dm,dim);
   DMPlexCreateFromDAG(dm, depth, numPoints, coneSize, cones, 
 coneOrientations, vertexCoords);

   DMPlexGetHeightStratum(dm, 0, , );
   DMPlexGetHeightStratum(dm, 1, , );

   for (int k =cStart;k>>>     double vol;
     double centroid[3];
     double normal[3];
     ierr = DMPlexComputeCellGeometryFVM(dm, k, , 
 centroid,NULL);CHKERRQ(ierr);
     printf("FVM: V=%f c=(%f %f %f) n=(%f %f 
 %f)\n",vol,centroid[0],centroid[1],centroid[2],
   normal[0],normal[1],normal[2]);
   }
   ierr = PetscFinalize();
   return 0;
 }



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



log3
Description: Binary data

Re: [petsc-users] Problems about SNES

2019-03-12 Thread Matthew Knepley via petsc-users
On Wed, Jan 16, 2019 at 10:59 PM Yingjie Wu via petsc-users <
petsc-users@mcs.anl.gov> wrote:

> Dear PETSc developers:
> Hi,
> During the process of testing the program, I found some questions about
> SNES. These are some basic questions that I have overlooked. Please help me
> to answer them.
> 1. Because my program uses - snes_mf_operator, there is no Jacobian
> matrix. Linear and non-linear step residuals are different in petsc. The
> linear step residuals are r_linear = J*δx-f(x). Since I don't have a
> Jacobian matrix, I don't know how to calculate the relative residuals of
> linear steps provided in petsc. Do we use the finite difference
> approximation matrix vector product when calculating the residuals?
>

PETSc is using a FD approximation to the Jacobian.


> 2. Read the user's manual for a brief introduction to the inexact Newton
> method, but I am very interested in the use of this method. I want to know
> how to use this method in petsc.
>

Inexact Newton is any method that does not solve the Newton system exactly,
so essentially any iterative solver. Do you mean quasi-Newton, which is
-snes_type qn


> 3. The default line search used by SNES in PETSc is bt, which often fails
> in program debugging. I don't know much about linesearch, and I'm curious
> to know why it failed. How can I supplement this knowledge?
>

Usually bt only fails if your search direction is not actually correct. Are
you sure your Jacobian is right?

  Thanks,

Matt


> Thanks,
> Yingjie
>


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


Re: [petsc-users] DMPlexComputeCellGeometryFVM: "Cannot handle faces with 1 vertices"

2019-03-12 Thread Chris Finn via petsc-users
Thank you,
I added DMPlexInterpolate(). Now I get a different error in DMPlexInterpolate():
...
[0]PETSC ERROR: Dimension 0 not supported
...

The source code and full output are attached. Anybody able to fix this?
regards
Chris

-- 
 Securely sent with Tutanota. Get your own encrypted, ad-free mailbox: 
 https://tutanota.com


Mar 12, 2019, 1:08 PM by knep...@gmail.com:

> On Tue, Mar 12, 2019 at 9:03 AM Chris Finn via petsc-users <> 
> petsc-users@mcs.anl.gov > > wrote:
>
>> Hello,
>> with the code below, I create a tetrahedron using DMPlexCreateFromDAG, then 
>> I try to run DMPlexComputeCellGeometryFVM on this cell. The latter call 
>> fails with output:
>>
>
> All the geometry stuff requires that you interpolate the mesh I think, Just 
> use DMPlexInterpolate().
>
>   Thanks,
>
>     Matt
>  
>
>>
>> ...
>> [0]PETSC ERROR: - Error Message -
>> [0]PETSC ERROR: Argument out of range
>> [0]PETSC ERROR: Cannot handle faces with 1 vertices
>> ...
>> (full output is attached).
>>
>> What am I doing wrong?
>> regards
>> Chris
>>
>> Here is the code:
>> static char help[] = "No help \n";
>>
>> #include 
>>
>> #undef __FUNCT__
>> #define __FUNCT__ "main"
>> int main(int argc,char **args)
>> {
>>   PetscErrorCode ierr;
>>   PetscInitialize(,,(char*)0,help);
>>
>>   DM    dm;
>>   int cStart,cEnd;
>>   int fStart,fEnd;
>>
>>   int depth = 3;
>>   int dim = 3;
>>
>>   PetscInt    numPoints[4]    = {1,4,6,4};
>>   PetscInt    coneSize[15] = {4,3,3,3,3,2,2,2,2,2,2,0,0,0,0};
>>   PetscInt    cones[28]    = {1,2,3,4, 5,9,8, 9,6,10, 10,8,7, 5,6,7, 
>> 11,12, 12,13, 13,11, 11,14, 12,14, 13,14};
>>   PetscInt    coneOrientations[28] = {0 };
>>   PetscScalar vertexCoords[12] = {0,0,0, 1,0,0, 0,1,0, 0,0,1};
>>
>>   DMCreate(PETSC_COMM_WORLD, );
>>   DMSetType(dm, DMPLEX);
>>   DMSetDimension(dm,dim);
>>   DMPlexCreateFromDAG(dm, depth, numPoints, coneSize, cones, 
>> coneOrientations, vertexCoords);
>>
>>   DMPlexGetHeightStratum(dm, 0, , );
>>   DMPlexGetHeightStratum(dm, 1, , );
>>
>>   for (int k =cStart;k>     double vol;
>>     double centroid[3];
>>     double normal[3];
>>     ierr = DMPlexComputeCellGeometryFVM(dm, k, , 
>> centroid,NULL);CHKERRQ(ierr);
>>     printf("FVM: V=%f c=(%f %f %f) n=(%f %f 
>> %f)\n",vol,centroid[0],centroid[1],centroid[2],
>>   normal[0],normal[1],normal[2]);
>>   }
>>   ierr = PetscFinalize();
>>   return 0;
>> }
>>
>>
>>
>
>
> -- 
> 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/ 
>

static char help[] = "No help \n";

#include 

#undef __FUNCT__
#define __FUNCT__ "main"
int main(int argc,char **args)
{
  PetscErrorCode ierr;
  PetscInitialize(,,(char*)0,help);

  DMdm;
  int cStart,cEnd;
  int fStart,fEnd;

  int depth = 3;
  int dim = 3;

  PetscIntnumPoints[4]= {1,4,6,4};
  PetscIntconeSize[15] = {4,3,3,3,3,2,2,2,2,2,2,0,0,0,0};
  PetscIntcones[28]= {1,2,3,4, 5,9,8, 9,6,10, 10,8,7, 5,6,7, 11,12, 12,13, 13,11, 11,14, 12,14, 13,14};
  PetscIntconeOrientations[28] = {0 };
  PetscScalar vertexCoords[12] = {0,0,0, 1,0,0, 0,1,0, 0,0,1};

  DMCreate(PETSC_COMM_WORLD, );
  DMSetType(dm, DMPLEX);
  DMSetDimension(dm,dim);
  ierr = DMPlexCreateFromDAG(dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords); CHKERRQ(ierr);


  DM idm=NULL; 
  ierr = DMPlexInterpolate(dm, );CHKERRQ(ierr);
  ierr = DMPlexCopyCoordinates(dm, idm);CHKERRQ(ierr);

  dm=idm;
  DMPlexGetHeightStratum(dm, 0, , );
  DMPlexGetHeightStratum(dm, 1, , );

  for (int k =cStart;k

log2
Description: Binary data


[petsc-users] DMPlexComputeCellGeometryFVM: "Cannot handle faces with 1 vertices"

2019-03-12 Thread Chris Finn via petsc-users
Hello,
with the code below, I create a tetrahedron using DMPlexCreateFromDAG, then I 
try to run DMPlexComputeCellGeometryFVM on this cell. The latter call fails 
with output:
...
[0]PETSC ERROR: - Error Message -
[0]PETSC ERROR: Argument out of range
[0]PETSC ERROR: Cannot handle faces with 1 vertices...
(full output is attached).

What am I doing wrong?
regards
Chris

Here is the code:
static char help[] = "No help \n";

#include 

#undef __FUNCT__
#define __FUNCT__ "main"
int main(int argc,char **args)
{
  PetscErrorCode ierr;
  PetscInitialize(,,(char*)0,help);

  DM    dm;
  int cStart,cEnd;
  int fStart,fEnd;

  int depth = 3;
  int dim = 3;

  PetscInt    numPoints[4]    = {1,4,6,4};
  PetscInt    coneSize[15] = {4,3,3,3,3,2,2,2,2,2,2,0,0,0,0};
  PetscInt    cones[28]    = {1,2,3,4, 5,9,8, 9,6,10, 10,8,7, 5,6,7, 
11,12, 12,13, 13,11, 11,14, 12,14, 13,14};
  PetscInt    coneOrientations[28] = {0 };
  PetscScalar vertexCoords[12] = {0,0,0, 1,0,0, 0,1,0, 0,0,1};

  DMCreate(PETSC_COMM_WORLD, );
  DMSetType(dm, DMPLEX);
  DMSetDimension(dm,dim);
  DMPlexCreateFromDAG(dm, depth, numPoints, coneSize, cones, coneOrientations, 
vertexCoords);

  DMPlexGetHeightStratum(dm, 0, , );
  DMPlexGetHeightStratum(dm, 1, , );

  for (int k =cStart;k

log
Description: Binary data


Re: [petsc-users] MatCompositeMerge + MatCreateRedundantMatrix

2019-03-12 Thread Marius Buerkle via petsc-users
I tried to follow your suggestions but it seems there is no MPICreateSubMatricesMPI for Fortran. Is this correct?

 



On Wed, Feb 20, 2019 at 6:57 PM Marius Buerkle  wrote:





ok, I think I understand now. I will give it a try and if there is some trouble comeback to you. thanks.




 

Cool.

 

   Matt

 




 

marius



 



On Tue, Feb 19, 2019 at 8:42 PM Marius Buerkle  wrote:




ok, so it seems there is no straight forward way to transfer data between PETSc matrices on different subcomms. Probably doing it by "hand" extracting the matricies on the subcomms create a MPI_INTERCOMM transfering the data to PETSC_COMM_WORLD and assembling them in a new PETSc matrix would be possible, right?



 

That sounds too complicated. Why not just reverse MPICreateSubMatricesMPI()? Meaning make it collective on the whole big communicator, so that you can swap out all the subcommunicator for the aggregation call, just like we do in that function.

Then its really just a matter of reversing the communication call.

 

   Matt 






 



On Tue, Feb 19, 2019 at 7:12 PM Marius Buerkle  wrote:





I see. This would work if the matrices are on different subcommumicators ? Is it possible to add this functionality ?




 

Hmm, no. That is specialized to serial matrices. You need the inverse of MatCreateSubMatricesMPI().

 

  Thanks,

 

     Matt

  




marius

 

 


You basically need the inverse of MatCreateSubmatrices(). I do not think we have that right now, but it could probably be done without too much trouble by looking at that code.
 

  Thanks,

 

     Matt

 


On Tue, Feb 19, 2019 at 6:15 AM Marius Buerkle via petsc-users  wrote:




Hi !

 

Is there some way to combine MatCompositeMerge with MatCreateRedundantMatrix? I basically want to create copies of a matrix from PETSC_COMM_WORLD to subcommunicators, do some work on each subcommunicator and than gather the results back to PETSC_COMM_WORLD, namely  I want to sum the  the invidual matrices from the subcommunicatos component wise and get the resulting matrix on PETSC_COMM_WORLD. Is this somehow possible without going through all the hassle of using MPI directly? 

 

marius




 

 
--







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/













 

 
--







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/














 

 
--







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/














 

 
--







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/












[petsc-users] Matrix Partitioning using PARMETIS

2019-03-12 Thread Eda Oktay via petsc-users
Hello,

I have a Laplacian matrix PL of matrix A and I try to partition A using
PARMETIS. Since PL is sequential and not adjacency matrix, I converted PL
to AL, then write the following code:

  ierr = MatConvert(PL,MATMPIADJ,MAT_INITIAL_MATRIX,);CHKERRQ(ierr);

  ierr = MatMeshToCellGraph(AL,2,);CHKERRQ(ierr);
  ierr = MatPartitioningCreate(MPI_COMM_WORLD,);CHKERRQ(ierr);
  ierr = MatPartitioningSetAdjacency(part,dual);CHKERRQ(ierr);

  ierr = MatPartitioningSetFromOptions(part);CHKERRQ(ierr);
  ierr = MatPartitioningApply(part,);CHKERRQ(ierr);
  ierr = ISView(partitioning,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);

  ierr = ISDestroy();CHKERRQ(ierr);
  ierr = MatPartitioningDestroy();CHKERRQ(ierr);

However, when I look at partitioning with ISView, the index set consists of
zeros only. Is that because I have only one processor and my codes are
written for only one processor, or is there another problem? I ran my code
with -mat_partitioning_type parmetis.

Thanks,

Eda