[petsc-dev] EM solver in XGC
Mark Adams
mfadams at lbl.gov
Tue Apr 14 13:53:38 CDT 2015
It looks like MatConvert is looking for "MatConvert_nest_aij_C"
but I see:
> git grep MatConvert_Nest
impls/nest/matnest.c:#define __FUNCT__ "MatConvert_Nest_AIJ"
impls/nest/matnest.c:PETSC_EXTERN PetscErrorCode MatConvert_Nest_AIJ(Mat
A,MatType newtype,MatReuse reuse,Mat *newmat)
Is there a problems with "Nest" vs "nest" here?
Mark
On Tue, Apr 14, 2015 at 2:46 PM, Mark Adams <mfadams at lbl.gov> wrote:
> I have tried
>
> call
> MatConvert(a_ts%FJacobian,MATMPIAIJ,MAT_INITIAL_MATRIX,a_ts%FJacobian2,ierr);CHKERRQ(ierr)
>
> and
>
> call
> MatConvert(a_ts%FJacobian,MATAIJ,MAT_INITIAL_MATRIX,a_ts%FJacobian2,ierr);CHKERRQ(ierr)
>
> in 'next' and it is not finding the specialized routine.
>
> I will look into it.
>
> Mark
>
>
> On Tue, Apr 14, 2015 at 1:41 PM, Matthew Knepley <knepley at gmail.com>
> wrote:
>
>> On Tue, Apr 14, 2015 at 12:40 PM, Matthew Knepley <knepley at gmail.com>
>> wrote:
>>
>>> On Tue, Apr 14, 2015 at 12:33 PM, Mark Adams <mfadams at lbl.gov> wrote:
>>>
>>>> MatNest does not support convert because it does not have a getrow.
>>>>
>>>
>>> I think this got fixed in dev. Anyhow, this supports my contention that
>>> MatNest is a piece of shit.
>>>
>>
>> Use 'next', it has
>>
>> PETSC_EXTERN PetscErrorCode MatConvert_Nest_AIJ(Mat A,MatType
>> newtype,MatReuse reuse,Mat *newmat)
>>
>> so it should work.
>>
>> Matt
>>
>>
>>>
>>>
>>>> Should we use a 4 block fieldsplit?
>>>>
>>>
>>> We give up our ability to check using exact solves, but go ahead and try
>>> it if you think it will make
>>> the physicists happy for now.
>>>
>>> Thanks,
>>>
>>> Matt
>>>
>>>
>>>> On Tue, Apr 14, 2015 at 1:26 PM, Mark Adams <mfadams at lbl.gov> wrote:
>>>>
>>>>> OK,
>>>>>
>>>>> On Tue, Apr 14, 2015 at 1:25 PM, Matthew Knepley <knepley at gmail.com>
>>>>> wrote:
>>>>>
>>>>>> On Tue, Apr 14, 2015 at 12:23 PM, Mark Adams <mfadams at lbl.gov> wrote:
>>>>>>
>>>>>>> I moved this to after the setup, which it has to do because it does
>>>>>>> not have the prefix or set options.
>>>>>>>
>>>>>>> It is not failing with not finding the IS. Will MatNest be able to
>>>>>>> merge two matrices together like this?
>>>>>>>
>>>>>>> I create a 4x4 MatNext but have just two fields in fieldSplit.
>>>>>>>
>>>>>>
>>>>>> Yep, I believe that will work. If not, just convert the MatNest to an
>>>>>> AIJ for now with MatConvert().
>>>>>>
>>>>>> Matt
>>>>>>
>>>>>>
>>>>>>> On Tue, Apr 14, 2015 at 1:18 PM, Matthew Knepley <knepley at gmail.com>
>>>>>>> wrote:
>>>>>>>
>>>>>>>> On Tue, Apr 14, 2015 at 12:13 PM, Mark Adams <mfadams at lbl.gov>
>>>>>>>> wrote:
>>>>>>>>
>>>>>>>>> I get an error where it looks like the number fields in field
>>>>>>>>> split is < 2. Do I need to set the size?
>>>>>>>>>
>>>>>>>>
>>>>>>>> No, each call to SetIS() just appends to the field list.
>>>>>>>>
>>>>>>>>
>>>>>>>>> I use this code but I do this before I have setup the TS. Should
>>>>>>>>> I do this after I setup the TS?
>>>>>>>>>
>>>>>>>>
>>>>>>>> This is why I hate using the API since it is difficult to get all
>>>>>>>> the setup done at the right time. First,
>>>>>>>> when you call PCFieldSplitSetIS(), is the PC type set?
>>>>>>>>
>>>>>>>> Thanks,
>>>>>>>>
>>>>>>>> Matt
>>>>>>>>
>>>>>>>>
>>>>>>>>> ! setup field split
>>>>>>>>> call TSGetSNES(a_ts,snes,ierr);CHKERRQ(ierr)
>>>>>>>>> call SNESGetKSP(snes,ksp,ierr);CHKERRQ(ierr)
>>>>>>>>> call KSPGetPC(ksp,pc,ierr);CHKERRQ(ierr)
>>>>>>>>>
>>>>>>>>> call
>>>>>>>>> MatGetSize(a_ts%mass,PETSC_NULL_INTEGER,npetsc,ierr);CHKERRQ(ierr)
>>>>>>>>> npetsc = 2*npetsc ! two blocks in each partition
>>>>>>>>> j = 0
>>>>>>>>> call ISCreateStride(w_comm,npetsc,j,ione,is,ierr);CHKERRQ(ierr)
>>>>>>>>> call PCFieldSplitSetIS(pc,'dt',is,ierr);CHKERRQ(ierr)
>>>>>>>>> call ISDestroy(is,ierr);CHKERRQ(ierr)
>>>>>>>>> j = npetsc
>>>>>>>>> call ISCreateStride(w_comm,npetsc,j,ione,is,ierr);CHKERRQ(ierr)
>>>>>>>>> call PCFieldSplitSetIS(pc,'aux',is,ierr);CHKERRQ(ierr)
>>>>>>>>> call ISDestroy(is,ierr);CHKERRQ(ierr)
>>>>>>>>>
>>>>>>>>>
>>>>>>>>> On Tue, Apr 14, 2015 at 12:24 PM, Matthew Knepley <
>>>>>>>>> knepley at gmail.com> wrote:
>>>>>>>>>
>>>>>>>>>> On Tue, Apr 14, 2015 at 11:16 AM, Mark Adams <mfadams at lbl.gov>
>>>>>>>>>> wrote:
>>>>>>>>>>
>>>>>>>>>>> Should I just use
>>>>>>>>>>>
>>>>>>>>>>> PetscErrorCode VecGetSubVector(Vec X,IS is,Vec *Y)
>>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>> Yep, that is exactly what I use it for.
>>>>>>>>>>
>>>>>>>>>> Matt
>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>>> ?
>>>>>>>>>>>
>>>>>>>>>>> On Tue, Apr 14, 2015 at 12:15 PM, Mark Adams <mfadams at lbl.gov>
>>>>>>>>>>> wrote:
>>>>>>>>>>>
>>>>>>>>>>>> I am thinking that I should create an IS that is my local part
>>>>>>>>>>>> of each field and then get a "sub vector" with that somehow. Looking into
>>>>>>>>>>>> this now.
>>>>>>>>>>>>
>>>>>>>>>>>> On Tue, Apr 14, 2015 at 12:00 PM, Mark Adams <mfadams at lbl.gov>
>>>>>>>>>>>> wrote:
>>>>>>>>>>>>
>>>>>>>>>>>>> Whoop, old code still there.
>>>>>>>>>>>>>
>>>>>>>>>>>>> I now have to get access to vectors somehow.
>>>>>>>>>>>>>
>>>>>>>>>>>>> I now get vectors for a field and do stuff with them like:
>>>>>>>>>>>>>
>>>>>>>>>>>>> subroutine FormIFunction(ts,t_dummy,a_XX,a_Xdot,a_FF,a_ts,ierr)
>>>>>>>>>>>>>
>>>>>>>>>>>>> call TSGetDM(ts,dm,ierr);CHKERRQ(ierr)
>>>>>>>>>>>>>
>>>>>>>>>>>>> call
>>>>>>>>>>>>> DMCompositeGetAccessArray(dm,a_FF,ifour,PETSC_NULL_INTEGER,Fsub,ierr);CHKERRQ(ierr)
>>>>>>>>>>>>> call
>>>>>>>>>>>>> DMCompositeGetAccessArray(dm,a_Xdot,ifour,PETSC_NULL_INTEGER,Xdotsub,ierr);CHKERRQ(ierr)
>>>>>>>>>>>>> one = 1d0
>>>>>>>>>>>>> call VecAXPY(Fsub(1),one,Xdotsub(1),ierr);CHKERRQ(ierr)
>>>>>>>>>>>>> call VecAXPY(Fsub(2),one,Xdotsub(2),ierr);CHKERRQ(ierr)
>>>>>>>>>>>>> call
>>>>>>>>>>>>> DMCompositeRestoreAccessArray(dm,a_Xdot,ifour,PETSC_NULL_INTEGER,Xdotsub,ierr);CHKERRQ(ierr)
>>>>>>>>>>>>> call
>>>>>>>>>>>>> DMCompositeRestoreAccessArray(dm,a_FF,ifour,PETSC_NULL_INTEGER,Fsub,ierr);CHKERRQ(ierr)
>>>>>>>>>>>>>
>>>>>>>>>>>>> call
>>>>>>>>>>>>> MatMultAdd(a_ts%FJacobian,a_XX,a_FF,a_FF,ierr);CHKERRQ(ierr)
>>>>>>>>>>>>>
>>>>>>>>>>>>> I need to get rid of DMCompositeGetAccessArray. What do I do?
>>>>>>>>>>>>>
>>>>>>>>>>>>> Thanks,
>>>>>>>>>>>>> Mark
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> On Tue, Apr 14, 2015 at 11:36 AM, Matthew Knepley <
>>>>>>>>>>>>> knepley at gmail.com> wrote:
>>>>>>>>>>>>>
>>>>>>>>>>>>>> On Tue, Apr 14, 2015 at 10:35 AM, Mark Adams <mfadams at lbl.gov
>>>>>>>>>>>>>> > wrote:
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>> And is this the solver (appended).
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> Do not call SNESSetDM(). We are not using DM at all in this
>>>>>>>>>>>>>> code.
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> Matt
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>> I am getting an error at this last line:
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>> .seealso: SNESGetDM(), KSPSetDM(), KSPGetDM()
>>>>>>>>>>>>>>> @*/
>>>>>>>>>>>>>>> PetscErrorCode SNESSetDM(SNES snes,DM dm)
>>>>>>>>>>>>>>> {
>>>>>>>>>>>>>>> PetscErrorCode ierr;
>>>>>>>>>>>>>>> KSP ksp;
>>>>>>>>>>>>>>> DMSNES sdm;
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>> PetscFunctionBegin;
>>>>>>>>>>>>>>> PetscValidHeaderSpecific(snes,SNES_CLASSID,1);
>>>>>>>>>>>>>>> if (dm) {ierr =
>>>>>>>>>>>>>>> PetscObjectReference((PetscObject)dm);CHKERRQ(ierr);}
>>>>>>>>>>>>>>> if (snes->dm) { /* Move the DMSNES context
>>>>>>>>>>>>>>> over to the new DM unless the new DM already has one */
>>>>>>>>>>>>>>> if (snes->dm->dmsnes && snes->dmAuto && !dm->dmsnes) {
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>> The problem is that the argument 'dm' to this method is NULL
>>>>>>>>>>>>>>> so I get a segv at the last term of this conditional.
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>> Mark
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>> -adv_ts_type beuler
>>>>>>>>>>>>>>> -adv_ts_monitor
>>>>>>>>>>>>>>> -adv_ts_view
>>>>>>>>>>>>>>> -adv_ts_max_step 1
>>>>>>>>>>>>>>> -adv_snes_max_it 1 # make >> 1 for
>>>>>>>>>>>>>>> nonlinear solve
>>>>>>>>>>>>>>> -adv_snes_lag_preconditioner -1 # make > 0 to update PC
>>>>>>>>>>>>>>> -adv_snes_lag_jacobian -1 # make > 0 for
>>>>>>>>>>>>>>> nonlinear Jacobian
>>>>>>>>>>>>>>> -adv_snes_monitor
>>>>>>>>>>>>>>> -adv_snes_converged_reason
>>>>>>>>>>>>>>> -adv_ksp_monitor
>>>>>>>>>>>>>>> -adv_ksp_converged_reason
>>>>>>>>>>>>>>> -adv_ksp_type richardson
>>>>>>>>>>>>>>> -adv_pc_type fieldsplit
>>>>>>>>>>>>>>> #-adv_pc_fieldsplit_type multiplicative
>>>>>>>>>>>>>>> -adv_pc_fieldsplit_type schur
>>>>>>>>>>>>>>> -adv_pc_fieldsplit_schur_fact_type full
>>>>>>>>>>>>>>> -adv_pc_fieldsplit_schur_precondition full
>>>>>>>>>>>>>>> #-adv_pc_fieldsplit_0_fields 0,1
>>>>>>>>>>>>>>> #-adv_pc_fieldsplit_1_fields 2,3
>>>>>>>>>>>>>>> -adv_fieldsplit_dt_pc_type lu
>>>>>>>>>>>>>>> -adv_fieldsplit_aux_pc_type lu
>>>>>>>>>>>>>>> -adv_fieldsplit_dt_pc_factor_mat_solver_package superlu
>>>>>>>>>>>>>>> -adv_fieldsplit_aux_pc_factor_mat_solver_package superlu
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> --
>>>>>>>>>>>>>> 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
>>>>>>>>>>
>>>>>>>>>
>>>>>>>>>
>>>>>>>>
>>>>>>>>
>>>>>>>> --
>>>>>>>> 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
>>>>>>
>>>>>
>>>>>
>>>>
>>>
>>>
>>> --
>>> 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-dev/attachments/20150414/319e1466/attachment.html>
More information about the petsc-dev
mailing list