On May 1, 2012, at 10:16 PM, Dmitry Karpeev wrote:

> Yeah, it's a bit of a hack right now -- DMCreateMatrix simply increfs and 
> return the matrix from the underlying NonlinearImplicitSystem.  The reason is 
> that initially the matrix is preallocated by unassembled, and MatDuplicate 
> doesn't want to duplicate an unassembled matrix.  I can fix that, but I'm a 
> bit puzzled by the fact that KSP messes 
> with the matrix.  Can you run with -snes_view and send me the output?  
> 
Dmitry,

Here it is:

------------- with the flag -----------------------

Linear solve converged due to CONVERGED_RTOL iterations 44
Calculating surface energy... done
Calculating fracture energy... done
elastic energy:13.74054,fracture energy:-0.00000,total:13.74054
applying variable bounds ... done
  NL step  0, |residual|_2 = 2.41046e-01
SNES Jacobian assembly... done
  Linear solve converged due to CONVERGED_RTOL iterations 3
  NL step  1, |residual|_2 = 5.10144e-07
  Linear solve did not converge due to DIVERGED_BREAKDOWN iterations 1
SNES Object: 1 MPI processes
  type: virs
  maximum iterations=50, maximum function evaluations=10000
  tolerances: relative=1e-08, absolute=1e-35, solution=1e-08
  total number of linear solver iterations=3
  total number of function evaluations=2
  KSP Object:   1 MPI processes
    type: gmres
      GMRES: restart=30, using Classical (unmodified) Gram-Schmidt 
Orthogonalization with no iterative refinement
      GMRES: happy breakdown tolerance 1e-30
    maximum iterations=10000, initial guess is zero
    tolerances:  relative=1e-05, absolute=1e-50, divergence=10000
    left preconditioning
    using PRECONDITIONED norm type for convergence test
  PC Object:   1 MPI processes
    type: ilu
      ILU: out-of-place factorization
      0 levels of fill
      tolerance for zero pivot 2.22045e-14
      using diagonal shift to prevent zero pivot
      matrix ordering: natural
      factor fill ratio given 1, needed 1
        Factored matrix follows:
          Matrix Object:           1 MPI processes
            type: seqaij
            rows=153, cols=153
            package used to perform factorization: petsc
            total: nonzeros=2401, allocated nonzeros=2401
            total number of mallocs used during MatSetValues calls =0
              not using I-node routines
    linear system matrix = precond matrix:
    Matrix Object:     1 MPI processes
      type: seqaij
      rows=153, cols=153
      total: nonzeros=2401, allocated nonzeros=2401
      total number of mallocs used during MatSetValues calls =0
        not using I-node routines
  SNESLineSearch Object:   1 MPI processes
    type: basic
    maxstep=1e+08, minlambda=1e-12
    tolerances: relative=1e-08, absolute=1e-15, lambda=1e-08
    maximum iterations=1
Fracture system solved at nonlinear iteration 1 , final nonlinear residual 
norm: 5.10144e-07

--- without the flag ------

applying variable bounds ... done
  NL step  0, |residual|_2 = 2.41046e-01
SNES Jacobian assembly... done
  Linear solve converged due to CONVERGED_RTOL iterations 3
  NL step  1, |residual|_2 = 5.10144e-07
SNES Jacobian assembly... done
  Linear solve converged due to CONVERGED_RTOL iterations 3
  NL step  2, |residual|_2 = 1.61153e-12
SNES Object: 1 MPI processes
  type: virs
  maximum iterations=50, maximum function evaluations=10000
  tolerances: relative=1e-08, absolute=1e-35, solution=1e-08
  total number of linear solver iterations=6
  total number of function evaluations=3
  KSP Object:   1 MPI processes
    type: gmres
      GMRES: restart=30, using Classical (unmodified) Gram-Schmidt 
Orthogonalization with no iterative refinement
      GMRES: happy breakdown tolerance 1e-30
    maximum iterations=10000, initial guess is zero
    tolerances:  relative=1e-05, absolute=1e-50, divergence=10000
    left preconditioning
    using PRECONDITIONED norm type for convergence test
  PC Object:   1 MPI processes
    type: ilu
      ILU: out-of-place factorization
      0 levels of fill
      tolerance for zero pivot 2.22045e-14
      using diagonal shift to prevent zero pivot
      matrix ordering: natural
      factor fill ratio given 1, needed 1
        Factored matrix follows:
          Matrix Object:           1 MPI processes
            type: seqaij
            rows=153, cols=153
            package used to perform factorization: petsc
            total: nonzeros=2401, allocated nonzeros=2401
            total number of mallocs used during MatSetValues calls =0
              not using I-node routines
    linear system matrix = precond matrix:
    Matrix Object:     1 MPI processes
      type: seqaij
      rows=153, cols=153
      total: nonzeros=2401, allocated nonzeros=2401
      total number of mallocs used during MatSetValues calls =0
        not using I-node routines
  SNESLineSearch Object:   1 MPI processes
    type: basic
    maxstep=1e+08, minlambda=1e-12
    tolerances: relative=1e-08, absolute=1e-15, lambda=1e-08
    maximum iterations=1
Fracture system solved at nonlinear iteration 2 , final nonlinear residual 
norm: 1.61153e-12



> Incidentally, there _was_ a reference counting bug that led to a breakdown 
> very much like what you reported.
> I'm attaching that patch and a small improvement to DMVariableBounds.  I need 
> to figure out how to go about creating matrices from DM_libMesh.
> 
 Let me know if you found anything or if I can be of any help.

Thanks,
Ata

> Dmitry.
> 
> On Tue, May 1, 2012 at 6:21 PM, Ataollah Mesgarnejad <[email protected]> 
> wrote:
> Dmitry,
> 
> Is there a way so SNES doesn't go calculate the Jacobian for every iteration 
> (my Jac does't change)? I used to use a flag to check if it's assembled but 
> when I do that with petsc-dm, KSP associated with SNES gives 
> DIVERGED_BREAKDOWN for the second iteration! looks like petsc-dm doesn't make 
> a copy of the matrix so when it's used in KSP solve it gets jumbled up.
> 
> Here is how it looks 
> 
> -with the flag:
> 
> applying variable bounds ... done
>   NL step  0, |residual|_2 = 8.90674e-04
> SNES Jacobian assembly... done
>   Linear solve converged due to CONVERGED_RTOL iterations 6
>   NL step  1, |residual|_2 = 2.90987e-05
>   Linear solve did not converge due to DIVERGED_BREAKDOWN iterations 1
> 
> this happens for every SNES solve.
> -----------------------
> 
> -without it:
> 
> applying variable bounds ... done
>   NL step  0, |residual|_2 = 6.08030e-02
> SNES hessian assembly... done
>   Linear solve converged due to CONVERGED_RTOL iterations 6
>   NL step  1, |residual|_2 = 2.61557e-07
> SNES hessian assembly... done
>   Linear solve converged due to CONVERGED_RTOL iterations 6
>   NL step  2, |residual|_2 = 1.35354e-12
> Fracture system solved at nonlinear iteration 2 , final nonlinear residual 
> norm: 1.35354e-12
> 
> ps: these are not from the same time step but you get the idea.
> 
> Best,
> Ata
> 
> On May 1, 2012, at 8:49 AM, Dmitry Karpeev wrote:
> 
>> As an fyi: I just built and ran miscellaneous_ex7 in parallel with the 
>> latest petsc-dev: 
>> http://petsc.cs.iit.edu/petsc/petsc-dev/rev/12d4a1ea0aca
>> If you decide you want to upgrade to this (or later) petsc-dev, you'd need a 
>> small patch for 
>> petsc_dm_nonlinear_solver.C to account for the latest API changes 
>> (attached).  
>> I am working on a bigger patch that adds more functionality to DM_libMesh 
>> and streamline
>> miscellaneous_ex7.  That will take some time to percolate up into the 
>> libMesh tree.  In the 
>> meantime you should be able to run with the small fix I've attached.
>> 
>> Let me know how it works.
>> Dmitry.
>> 
>> 
>> On Tue, May 1, 2012 at 8:35 AM, Ataollah Mesgarnejad 
>> <[email protected]> wrote:
>> Dmitry,
>> 
>> I believe that the error is indeed due to petsc-dev version. I'm going to 
>> check this during the next few hours and get back to you.
>> 
>> Thanks,
>> Ata
>> 
>> 
>> On Tue, May 1, 2012 at 8:29 AM, Dmitry Karpeev <[email protected]> wrote:
>> Ata,
>> I can't reproduce the problem: miscellaneous_ex7 runs fine for me on 2 procs 
>> with a relatively recent version 
>> of petsc-dev.  Can you send me the full PETSc error output?
>> 
>> Thanks.
>> Dmitry.
>> 
>> 
>> On Mon, Apr 30, 2012 at 8:56 PM, Ataollah Mesgarnejad 
>> <[email protected]> wrote:
>> 
>> On Apr 30, 2012, at 8:39 PM, Dmitry Karpeev wrote:
>> 
>>> Ata,
>>> example7 segfaults when you run in parallel?
>>> Does it do it on 2 procs?
>> Yes and yes. 
>> 
>> here is the stack:
>> (gdb) where
>> #0  0x00000038fbaaada0 in __nanosleep_nocancel () from /lib64/libc.so.6
>> #1  0x00000038fbaaac30 in sleep () from /lib64/libc.so.6
>> #2  0x00002ae457b35c40 in PetscSleep(double) ()
>>    from /share/apps/petsc-dev/Linux-intel12.1-mef90-O/lib/libpetsc.so
>> #3  0x00002ae457b7dffa in PetscAttachDebugger() ()
>>    from /share/apps/petsc-dev/Linux-intel12.1-mef90-O/lib/libpetsc.so
>> #4  0x00002ae457b7e95f in PetscAttachDebuggerErrorHandler(int, int, char 
>> const*, char const*, char const*, int, PetscErrorType, char const*, void*) ()
>>    from /share/apps/petsc-dev/Linux-intel12.1-mef90-O/lib/libpetsc.so
>> #5  0x00002ae457b817be in PetscError(int, int, char const*, char const*, 
>> char const*, int, PetscErrorType, char const*, ...) ()
>>    from /share/apps/petsc-dev/Linux-intel12.1-mef90-O/lib/libpetsc.so
>> #6  0x00002ae457b82ea7 in PetscDefaultSignalHandler(int, void*) ()
>>    from /share/apps/petsc-dev/Linux-intel12.1-mef90-O/lib/libpetsc.so
>> #7  0x00002ae457b82883 in PetscSignalHandler_Private ()
>>    from /share/apps/petsc-dev/Linux-intel12.1-mef90-O/lib/libpetsc.so
>> #8  <signal handler called>
>> #9  0x00002ae458189340 in DMCoarsenHookAdd(_p_DM*, int (*)(_p_DM*, _p_DM*, 
>> void*), int (*)(_p_DM*, _p_Mat*, _p_Vec*, _p_Mat*, _p_DM*, void*), void*) ()
>>    from /share/apps/petsc-dev/Linux-intel12.1-mef90-O/lib/libpetsc.so
>> #10 0x00002ae45835f431 in DMSNESGetContext(_p_DM*, _n_SNESDM**) ()
>>    from /share/apps/petsc-dev/Linux-intel12.1-mef90-O/lib/libpetsc.so
>> #11 0x00002ae45835ea34 in DMSNESGetFunction(_p_DM*, int (**)(_p_SNES*, 
>> _p_Vec*, ---Type <return> to continue, or q <return> to quit---
>> _p_Vec*, void*), void**) ()
>>    from /share/apps/petsc-dev/Linux-intel12.1-mef90-O/lib/libpetsc.so
>> #12 0x00002ae458393d1c in SNESGetFunction(_p_SNES*, _p_Vec**, int 
>> (**)(_p_SNES*, _p_Vec*, _p_Vec*, void*), void**) ()
>>    from /share/apps/petsc-dev/Linux-intel12.1-mef90-O/lib/libpetsc.so
>> #13 0x00002ae45839b9eb in SNESSetUp(_p_SNES*) ()
>>    from /share/apps/petsc-dev/Linux-intel12.1-mef90-O/lib/libpetsc.so
>> #14 0x00002ae45839a6b0 in SNESSolve(_p_SNES*, _p_Vec*, _p_Vec*) ()
>>    from /share/apps/petsc-dev/Linux-intel12.1-mef90-O/lib/libpetsc.so
>> 
>> #15 0x00000000010976ab in 
>> libMesh::PetscDMNonlinearSolver<libMesh::Number>::solve (this=0x35506d0, 
>> jac_in=..., x_in=..., r_in=...)
>>     at src/solvers/petsc_dm_nonlinear_solver.C:450
>> #16 0x0000000000cd604f in libMesh::NonlinearImplicitSystem::solve (
>>     this=0x353b7e0) at src/systems/nonlinear_implicit_system.C:178
>> #17 0x0000000000426082 in Biharmonic::step (this=0x34f2000, 
>> dt=@0x7fff95edca28)
>>     at biharmonic.C:110
>> #18 0x0000000000426373 in Biharmonic::run (this=0x34f2000) at 
>> biharmonic.C:152
>> warning: Range for type (null) has invalid bounds 0..-111
>> #19 0x0000000000436cf5 in main (argc=12, argv=0x7fff95edcc18)
>>     at miscellaneous_ex7.C:64
>> 
>> Ata
>> 
>>> 
>>> Dmitry.
>>> 
>>> On Mon, Apr 30, 2012 at 8:31 PM, Ataollah Mesgarnejad 
>>> <[email protected]> wrote:
>>> Dmitry,
>>> 
>>> I had a very good experience on single processor but on multiple processors 
>>> I get segmentation error at line where we call solve (I tried both my code 
>>> and  example7) . I tried your example without -vi --use-petsc-dm and it 
>>> works just fine. Right now I'm trying to compile a debug version and look 
>>> into it; meanwhile I thought it would be a good idea if I let you know.
>>> 
>>> Thanks,
>>> Ata
>>> 
>>> On Apr 30, 2012, at 1:30 PM, Dmitry Karpeev wrote:
>>> 
>>>> 
>>>> 
>>>> On Mon, Apr 30, 2012 at 1:28 PM, Ataollah Mesgarnejad 
>>>> <[email protected]> wrote:
>>>> 
>>>> On Apr 30, 2012, at 1:13 PM, Dmitry Karpeev wrote:
>>>> 
>>>>> 
>>>>> 
>>>>> On Mon, Apr 30, 2012 at 1:03 PM, Ataollah Mesgarnejad 
>>>>> <[email protected]> wrote:
>>>>> Dmitry,
>>>>> 
>>>>> I had time to come back and work on the SNESVI/Libmesh program. The way 
>>>>> you setup your code is particularly hard for me to follow in my own code 
>>>>> since it requires a significant amount of  extra work (to define classes 
>>>>> and setup the inheritances..). What I did instead was to try add a 
>>>>> TransientNonlinearSystem to my EquationSystems object:
>>>>> 
>>>>>   TransientNonlinearImplicitSystem & fracture_system =            
>>>>> equation_system.add_system<TransientNonlinearImplicitSystem>("fracture_system");
>>>>> 
>>>>> and then to add the bounds and Objective and Residual functions:
>>>>> 
>>>>>   fracture_system.attach_init_function(init_ic_fracture);
>>>>>   fracture_system.nonlinear_solver->jacobian= PsiFormJacobian;
>>>>>   fracture_system.nonlinear_solver->residual= PsiFormResidual;
>>>>>   fracture_system.nonlinear_solver->bounds = PsiFormBounds;
>>>>> 
>>>>> what I see now is that SNES solve does not apply the bounds. So here is 
>>>>> my questions:
>>>>> 
>>>>> 1-  Is there inherently wrong with setting up the system specially the 
>>>>> bounds the easy way like I did?
>>>>> Ata,
>>>>> Doing things this "simple" way is just fine: there is no reason to follow 
>>>>> the relatively convoluted structure I used in that 
>>>>> example.  When I have time, I'll rewrite it to simplify things. 
>>>>> 2- Where is the right SNES type is chosen ? I guess this should be hard 
>>>>> wired somewhere in NonlinearImplicitSystem!
>>>>> 3- When I use -snes_view (with your example or my own code ) I get 
>>>>> snes_type ls? shouldn't this be VI ?
>>>>>  
>>>>> What's likely missing is the command-line option --use-petsc-dm.  It's 
>>>>> necessary since libMesh wants to preserve the usual behavior as the 
>>>>> default.  In particular, because SNESVI and the DM hookup aren't really 
>>>>> possible prior to petsc-3.2.
>>>>> In fact, it might not be possible without a recent petsc-dev (let me 
>>>>> verify that, once I get to my computer).
>>>>> 
>>>>> In the meantime, could you give --use-petsc-dm a try?
>>>>> 
>>>> 
>>>> I just saw that in the example help line and used it. It works like a 
>>>> charm.
>>>> Great!
>>>> Dmitry. 
>>>> I'll go on check for the convergence and see what I can get.
>>>> 
>>>> Thanks,
>>>> Ata
>>>>> Thanks.
>>>>> 
>>>>> Dmitry.
>>>>> 
>>>>> Thanks,
>>>>> Ata
>>>>>  
>>>>> On Apr 12, 2012, at 10:00 AM, Dmitry Karpeev wrote:
>>>>> 
>>>>>> Convergence is something we might need to work on -- I noticed some 
>>>>>> chatter in a simple bounded biharmonic heat equation (a stripped down 
>>>>>> version of Cahn-Hilliard).  
>>>>>> Thanks.
>>>>>> Dmitry.
>>>>>> 
>>>>>> On Thu, Apr 12, 2012 at 9:55 AM, Ataollah Mesgarnejad 
>>>>>> <[email protected]> wrote:
>>>>>> Thanks Dmitry,
>>>>>> 
>>>>>> It would be wonderful. I made a wrapper for SNESVI in myself for my 
>>>>>> program but had difficulty achieving good convergence. Right now I am 
>>>>>> occupied with other issues but I'll get back to this as soon as I get a 
>>>>>> chance.
>>>>>> 
>>>>>> Ata
>>>>>> 
>>>>>> On Apr 10, 2012, at 9:55 PM, Dmitry Karpeev wrote:
>>>>>> 
>>>>>>> Ata,
>>>>>>> There is now a way to use SNESVI from libMesh.
>>>>>>> If interested, pull a recent revision from sourceforge and take a look 
>>>>>>> at examples/miscellaneous/miscellaneous_ex7.
>>>>>>> Let me know, if you have any questions.
>>>>>>> 
>>>>>>> Thanks.
>>>>>>> Dmitry.
>>>>>>> 
>>>>>>> On Wed, Dec 14, 2011 at 12:45 PM, Dmitry Karpeev <[email protected]> 
>>>>>>> wrote:
>>>>>>> Okay, thanks for your feedback.
>>>>>>> I'll keep you updated on this work.
>>>>>>> 
>>>>>>> Dmitry.
>>>>>>> 
>>>>>>> 
>>>>>>> On Wed, Dec 14, 2011 at 12:42 PM, Ataollah Mesgarnejad 
>>>>>>> <[email protected]> wrote:
>>>>>>> 
>>>>>>> 
>>>>>>> On Wed, Dec 14, 2011 at 12:35 PM, Dmitry Karpeev <[email protected]> 
>>>>>>> wrote:
>>>>>>> 
>>>>>>> 
>>>>>>> On Wed, Dec 14, 2011 at 12:19 PM, Ataollah Mesgarnejad 
>>>>>>> <[email protected]> wrote:
>>>>>>> 
>>>>>>> 
>>>>>>> On Wed, Dec 14, 2011 at 11:59 AM, Dmitry Karpeev <[email protected]> 
>>>>>>> wrote:
>>>>>>> Ata,
>>>>>>> 
>>>>>>> A couple of things: 
>>>>>>> 1. SNESVI solves a VI under what's known as "box constraints": the Vec 
>>>>>>> x that the residual operates on is entrywise between Vecs xlow and 
>>>>>>> xhigh.  In order to enable 
>>>>>>> constraints on a linear combination of the Vec entries (as you 
>>>>>>> suggest), the problem
>>>>>>> has to be augmented with extra variables -- these linear combinations.  
>>>>>>> This, of course,
>>>>>>> changes the Jacobian, preconditioner, etc.  The user provides only the 
>>>>>>> Jacobian for the original system, so libMesh would have to 
>>>>>>> automatically construct the Jacobian for the extended system, and 
>>>>>>> precondition it using the user provided PC for the original subsystem. 
>>>>>>> This is not necessarily the unreasonable thing to do in many 
>>>>>>> circumstances, but will take more thinking than a straighforward nodal 
>>>>>>> constraint.  The user API for specifying the constraints, however, 
>>>>>>> becomes rather straighforward.
>>>>>>> 
>>>>>>> 2. In terms of FEM convergence, constraining at the nodes is sufficient 
>>>>>>> (see, for example, 
>>>>>>> Glowinski's book: 
>>>>>>> http://books.google.com/books/about/Numerical_analysis_of_variational_inequa.html?id=Pf4ed2mtbx4C)
>>>>>>> You can show that if you impose the nodal constraints, in the limit of 
>>>>>>> mesh size going to 0 the inequality is satisfied a.e.
>>>>>>> 
>>>>>>> In this case no funny manipulation of the Jacobian or the PC is 
>>>>>>> necessary, but the user API to specify the constraint is a bit 
>>>>>>> amiguous: how do you constrain the node?  From a containing element?  
>>>>>>> Then there are multiple, possibly incompatible, constraints.
>>>>>>> Or you can simply go through all the nodes in the system once using 
>>>>>>> node iterator. Here is how I think you can do it:
>>>>>>> 
>>>>>>> -------------------------------------
>>>>>>> //make two PetscVector ub,lb from two variables you already added to 
>>>>>>> the system
>>>>>>> NumericVector<Number> &ub_sys = fracture_system.get_vector("ub");
>>>>>>> NumericVector<Number> &lb_sys = fracture_system.get_vector("lb");
>>>>>>> PetscVector<Number>&lb  =  
>>>>>>> libmesh_cast_ptr<PetscVector<Number>*>(&ub_sys);
>>>>>>> PetscVector<Number>&ub  =  
>>>>>>> libmesh_cast_ptr<PetscVector<Number>*>(&lb_sys);
>>>>>>> 
>>>>>>> then in the bound routine do:
>>>>>>> 
>>>>>>> // Get a constant reference to the mesh object.
>>>>>>> const MeshBase& mesh = es.get_mesh();
>>>>>>> //Get the iterator
>>>>>>> MeshBase::const_node_iterator       nd     = mesh.local_nodes_begin();
>>>>>>> const MeshBase::const_node_iterator end_nd = mesh.local_nodes_end();
>>>>>>>         for (; nd != end_nd; nd++){
>>>>>>>                 const Node* node = *nd;
>>>>>>>                 //get the degree of freedom of node 
>>>>>>>                unsigned int dn = node->dof_number(0,0,0); // system 0 , 
>>>>>>> var 0, component 0
>>>>>>>                //set the ub,lb
>>>>>>>                ub(dn)= A;
>>>>>>>                lb(dn) = B;
>>>>>>>  }
>>>>>>> //close the vectors
>>>>>>> ub.close();
>>>>>>> lb.close();
>>>>>>> 
>>>>>>> You could do that. That's probably the mode we can support almost 
>>>>>>> immediately.
>>>>>>> The issue here is that you lose some of the "context" by processing 
>>>>>>> each node by
>>>>>>> itself rather as part of  an element.  For example what if each element 
>>>>>>> is a different material?
>>>>>>> I guess in that case that property should be set on the nodes as well.
>>>>>>> 
>>>>>>> Exactly one can always pass mat prop as a Constant MONOMIAL (e.g. 
>>>>>>> system.add_variable("mat_prop",CONSTANT,MONOMIAL);)  variable with 
>>>>>>> values on each node.  
>>>>>>> 
>>>>>>> Ata
>>>>>>>  
>>>>>>> Dmitry.
>>>>>>> -----------------------------------------------------------------
>>>>>>> Ata
>>>>>>>  
>>>>>>> 
>>>>>>> Dmitry.
>>>>>>> 
>>>>>>> 
>>>>>>> On Wed, Dec 14, 2011 at 11:17 AM, Ataollah Mesgarnejad 
>>>>>>> <[email protected]> wrote:
>>>>>>> I hope I understood the whole thing correctly. 
>>>>>>> 
>>>>>>> On Wed, Dec 14, 2011 at 10:56 AM, Dmitry Karpeev <[email protected]> 
>>>>>>> wrote:
>>>>>>> 
>>>>>>> 
>>>>>>> On Wed, Dec 14, 2011 at 10:50 AM, Ataollah Mesgarnejad 
>>>>>>> <[email protected]> wrote:
>>>>>>> I don't think setting the variable bound on each node separately is a 
>>>>>>> good idea in general.
>>>>>>> I'm not sure what you mean by this.  Could you clarify?  I don't think 
>>>>>>> I'm advocating this. 
>>>>>>>   
>>>>>>> Bounds should be set as a sum on each elements degrees of freedom; this 
>>>>>>> way there should not be any nodal inconsistencies.
>>>>>>> I'm not sure what you mean here either.  A bound is usually not a sum 
>>>>>>> (e.g., a bound on a concentration is 0 <= c <=1).
>>>>>>> Perhaps your application expresses the total bound as a sum?  I'd be 
>>>>>>> very interested in learning more about that case.
>>>>>>>  
>>>>>>> OK here is my thought process: one optimally wants to bound the 
>>>>>>> variable over the whole Element domain between lb,ub. So we would have 
>>>>>>> something like this:
>>>>>>> 
>>>>>>> lb \leq \sum_{i=1}^{n=NDOF} \phi_{i}(x) \hat{u}_{i} \leq ub, \forall x 
>>>>>>> \in El                                                       (Eq-1)
>>>>>>> 
>>>>>>> so it's not exactly that the nodal values should be bounded 
>>>>>>> (\hat{u}_{i}) unless it's a FD application where \phi_i(x) is delta 
>>>>>>> function !! That said I think given a normalized basis, bounding nodal 
>>>>>>> values is stricter requirement that (Eq-1). if that's what you meant by 
>>>>>>> inconsistencies between node bounds then the stronger bound would be my 
>>>>>>> first priority. 
>>>>>>>  
>>>>>>> Ata
>>>>>>> 
>>>>>>>  
>>>>>>> Of course it might run into a situations with semi-smooth VI.
>>>>>>> The type of VI solver is independent of how the bounds are specified: 
>>>>>>> once we know the upper and lower bounds for each element of a Vec, we 
>>>>>>> can solve the VI with either a semismooth method or an active set 
>>>>>>> method.
>>>>>>> 
>>>>>>> I really appreciate your feedback.  I hope we can get to the point 
>>>>>>> where things are as clear as possible so that you can
>>>>>>> profitably use the new capabilities.
>>>>>>> 
>>>>>>> Dmitry. 
>>>>>>> 
>>>>>>> Thanks,
>>>>>>> Ata
>>>>>>> 
>>>>>>> On Wed, Dec 14, 2011 at 9:47 AM, Dmitry Karpeev <[email protected]> 
>>>>>>> wrote:
>>>>>>> Are you familiar with SNESVI?  The main thing I'm looking at right now 
>>>>>>> is the libMesh API for specifying
>>>>>>> the bounds on the variables. It seems to me that the user would ideally 
>>>>>>> want to specify the bounds on the variable
>>>>>>> one element at a time, since that's the default traversal mode for 
>>>>>>> libMesh.  Then, of course, many nodal values
>>>>>>> are traversed many time, allowing for an inconsistency.  The question 
>>>>>>> is, what you, as a user, expect in that situation?
>>>>>>> Should the latest  bound that the user has set be applied or the 
>>>>>>> tightest of the bounds the user has supplied?
>>>>>>> 
>>>>>>> Dmitry.
>>>>>>> 
>>>>>>> 
>>>>>>> On Wed, Dec 14, 2011 at 9:28 AM, Ataollah Mesgarnejad 
>>>>>>> <[email protected]> wrote:
>>>>>>> Thanks Dmitry, I look forward to use it. Also if it's available as a 
>>>>>>> dev version I would love to try to work with it.
>>>>>>> 
>>>>>>> Ata
>>>>>>> 
>>>>>>> 
>>>>>>> On Wed, Dec 14, 2011 at 8:24 AM, Dmitry Karpeev <[email protected]> 
>>>>>>> wrote:
>>>>>>> We (ANL) are working on an patch to libMesh to enable this 
>>>>>>> functionality.
>>>>>>> The motivation is to use it through MOOSE, but it will be exposed as a 
>>>>>>> libMesh API.
>>>>>>> Dmitry.
>>>>>>> 
>>>>>>> On Wed, Dec 14, 2011 at 7:53 AM, Ataollah Mesgarnejad 
>>>>>>> <[email protected]> wrote:
>>>>>>> Dear all,
>>>>>>> 
>>>>>>> I wonder if there is compatible versionIs it possible to
>>>>>>> use PETScNonlinearSolver class with SNESVI (PETSc/3.2)  to solve a
>>>>>>> variational inequality problem?
>>>>>>> 
>>>>>>> Thanks,
>>>>>>> Ata Mesgarnejad
>>>>>>> ------------------------------------------------------------------------------
>>>>>>> Cloud Computing - Latest Buzzword or a Glimpse of the Future?
>>>>>>> This paper surveys cloud computing today: What are the benefits?
>>>>>>> Why are businesses embracing it? What are its payoffs and pitfalls?
>>>>>>> http://www.accelacomm.com/jaw/sdnl/114/51425149/
>>>>>>> _______________________________________________
>>>>>>> Libmesh-users mailing list
>>>>>>> [email protected]
>>>>>>> https://lists.sourceforge.net/lists/listinfo/libmesh-users
>>>>>>> 
>>>>>>> 
>>>>>>> 
>>>>>>> 
>>>>>>> -- 
>>>>>>> A. Mesgarnejad
>>>>>>> PhD Student, Research Assistant
>>>>>>> Mechanical Engineering Department
>>>>>>> Louisiana State University
>>>>>>> 2203 Patrick F. Taylor Hall
>>>>>>> Baton Rouge, La 70803
>>>>>>> 
>>>>>>> 
>>>>>>> 
>>>>>>> 
>>>>>>> -- 
>>>>>>> A. Mesgarnejad
>>>>>>> PhD Student, Research Assistant
>>>>>>> Mechanical Engineering Department
>>>>>>> Louisiana State University
>>>>>>> 2203 Patrick F. Taylor Hall
>>>>>>> Baton Rouge, La 70803
>>>>>>> 
>>>>>>> 
>>>>>>> 
>>>>>>> 
>>>>>>> 
>>>>>>> 
>>>>>>> -- 
>>>>>>> A. Mesgarnejad
>>>>>>> PhD Student, Research Assistant
>>>>>>> Mechanical Engineering Department
>>>>>>> Louisiana State University
>>>>>>> 2203 Patrick F. Taylor Hall
>>>>>>> Baton Rouge, La 70803
>>>>>>> 
>>>>>>> 
>>>>>>> 
>>>>>>> 
>>>>>>> -- 
>>>>>>> A. Mesgarnejad
>>>>>>> PhD Student, Research Assistant
>>>>>>> Mechanical Engineering Department
>>>>>>> Louisiana State University
>>>>>>> 2203 Patrick F. Taylor Hall
>>>>>>> Baton Rouge, La 70803
>>>>>>> 
>>>>>>> 
>>>>>> 
>>>>>> 
>>>>> 
>>>>> 
>>>> 
>>>> 
>>> 
>>> 
>> 
>> 
>> 
>> 
>> 
>> -- 
>> A. Mesgarnejad
>> PhD Student, Research Assistant
>> Mechanical Engineering Department
>> Louisiana State University
>> 2203 Patrick F. Taylor Hall
>> Baton Rouge, La 70803
>> 
>> <petsc_dm.patch>
> 
> 
> <0001-PETSc-DM-incref-bug-fix.patch><0002-PETSc-DM-initializes-bounds-with-Inf-now.patch>

------------------------------------------------------------------------------
Live Security Virtual Conference
Exclusive live event will cover all the ways today's security and 
threat landscape has changed and how IT managers can respond. Discussions 
will include endpoint security, mobile security and the latest in malware 
threats. http://www.accelacomm.com/jaw/sfrnl04242012/114/50122263/
_______________________________________________
Libmesh-users mailing list
[email protected]
https://lists.sourceforge.net/lists/listinfo/libmesh-users

Reply via email to