On Tue, Apr 14, 2015 at 5:19 PM, Mark Adams <[email protected]> wrote: > Matt, in looking at MatConvert_Nest_AIJ, the code looks like it is 2 years > old. > > I need a MatConvert_Nest_MPIAIJ, should I write this or does it exist > someplace? >
This is MPIAIJ. Look at the code. Matt > > Mark > > > > On Tue, Apr 14, 2015 at 2:54 PM, Mark Adams <[email protected]> wrote: > >> Also, this is with MPI so I need MPIAIJ. Is that supposed to work? >> >> On Tue, Apr 14, 2015 at 2:53 PM, Mark Adams <[email protected]> wrote: >> >>> 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 <[email protected]> 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 <[email protected]> >>>> wrote: >>>> >>>>> On Tue, Apr 14, 2015 at 12:40 PM, Matthew Knepley <[email protected]> >>>>> wrote: >>>>> >>>>>> On Tue, Apr 14, 2015 at 12:33 PM, Mark Adams <[email protected]> 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 <[email protected]> wrote: >>>>>>> >>>>>>>> OK, >>>>>>>> >>>>>>>> On Tue, Apr 14, 2015 at 1:25 PM, Matthew Knepley <[email protected] >>>>>>>> > wrote: >>>>>>>> >>>>>>>>> On Tue, Apr 14, 2015 at 12:23 PM, Mark Adams <[email protected]> >>>>>>>>> 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 < >>>>>>>>>> [email protected]> wrote: >>>>>>>>>> >>>>>>>>>>> On Tue, Apr 14, 2015 at 12:13 PM, Mark Adams <[email protected]> >>>>>>>>>>> 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 < >>>>>>>>>>>> [email protected]> wrote: >>>>>>>>>>>> >>>>>>>>>>>>> On Tue, Apr 14, 2015 at 11:16 AM, Mark Adams <[email protected]> >>>>>>>>>>>>> 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 <[email protected] >>>>>>>>>>>>>> > 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 < >>>>>>>>>>>>>>> [email protected]> 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 < >>>>>>>>>>>>>>>> [email protected]> wrote: >>>>>>>>>>>>>>>> >>>>>>>>>>>>>>>>> On Tue, Apr 14, 2015 at 10:35 AM, Mark Adams < >>>>>>>>>>>>>>>>> [email protected]> 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 >>>>> >>>> >>>> >>> >> > -- 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
