On Wed, 2 Mar 2011, Boyce Griffith wrote:

> On 3/1/11 6:41 PM, Roy Stogner wrote:
>> 
>> On Thu, 2 Dec 2010, Boyce Griffith wrote:
>> 
>>> The patched code does the "right thing" only when
>>> asymmetric_constraint_rows=true when calling
>>> DofMap::constrain_element_{matrix,vector,matrix_and_vector}(). This is
>>> the default, but if asymmetric_constraint_rows=false, the code will
>>> silently do the "wrong thing," probably yielding bizarre results. One
>>> approach would be to keep track of whether a particular constraint row
>>> is provided by the user or not, and to require
>>> asymmetric_constraint_rows==true when such constraints are present.
>> 
>> Thanks for the heads-up. I've got some FEMSystem-based codes that
>> could probably be adjusted pretty easily to test that case, which we'd
>> want to fix before committing.
>
> Let me know if you want me to think about what to do for this case, or if it 
> isn't clear what goes wrong if asymmetric_constraint_rows=false.

I haven't looked at the problem too closely yet, but it might be
better to think about the math before the code.

For the inhomogeneous cases: we basically turn K*x = F, x=C*y into
C^T*K*C*y=C^T*F.  That's still a symmetric matrix.  It's a smaller
symmetric matrix in the case where length(y) < length(x), but we don't
actually condense the matrix, so we end up having indeterminate
hanging node entries in y.  In the symmetric case we make the matrix a
diagonal for those entries, and then we fix them up later after the
linear solve with enforce_constraints_exactly.  In the asymmetric case
we reuse those entries' lines with the constraint equations so that
the solver will enforce them (if only inexactly).

So now we've got x=C*y+H

Which means if we do things analogously in the symmetric case we want

C^T*K*(C*y+H)=C^T*F
C^T*K*C*y = C^T*(F-K*H)

Hmm... this isn't what your patch does, and when I think about it we
probably *don't* want to do things analogously (it's nice that in the
homogeneous case x is also just y with some entries omitted, but the
"x=C*y+H" model doesn't preserve that property).  Is there any
similarly clean way to formulate what your patch is doing?
---
Roy

------------------------------------------------------------------------------
Free Software Download: Index, Search & Analyze Logs and other IT data in 
Real-Time with Splunk. Collect, index and harness all the fast moving IT data 
generated by your applications, servers and devices whether physical, virtual
or in the cloud. Deliver compliance at lower cost and gain new business 
insights. http://p.sf.net/sfu/splunk-dev2dev 
_______________________________________________
Libmesh-devel mailing list
Libmesh-devel@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/libmesh-devel

Reply via email to