You can subscribe to this list here.
2003 
_{Jan}

_{Feb}

_{Mar}

_{Apr}

_{May}

_{Jun}

_{Jul}

_{Aug}

_{Sep}
(2) 
_{Oct}
(2) 
_{Nov}
(27) 
_{Dec}
(31) 

2004 
_{Jan}
(6) 
_{Feb}
(15) 
_{Mar}
(33) 
_{Apr}
(10) 
_{May}
(46) 
_{Jun}
(11) 
_{Jul}
(21) 
_{Aug}
(15) 
_{Sep}
(13) 
_{Oct}
(23) 
_{Nov}
(1) 
_{Dec}
(8) 
2005 
_{Jan}
(27) 
_{Feb}
(57) 
_{Mar}
(86) 
_{Apr}
(23) 
_{May}
(37) 
_{Jun}
(34) 
_{Jul}
(24) 
_{Aug}
(17) 
_{Sep}
(50) 
_{Oct}
(24) 
_{Nov}
(10) 
_{Dec}
(60) 
2006 
_{Jan}
(47) 
_{Feb}
(46) 
_{Mar}
(127) 
_{Apr}
(19) 
_{May}
(26) 
_{Jun}
(62) 
_{Jul}
(47) 
_{Aug}
(51) 
_{Sep}
(61) 
_{Oct}
(42) 
_{Nov}
(50) 
_{Dec}
(33) 
2007 
_{Jan}
(60) 
_{Feb}
(55) 
_{Mar}
(77) 
_{Apr}
(102) 
_{May}
(82) 
_{Jun}
(102) 
_{Jul}
(169) 
_{Aug}
(117) 
_{Sep}
(80) 
_{Oct}
(37) 
_{Nov}
(51) 
_{Dec}
(43) 
2008 
_{Jan}
(71) 
_{Feb}
(94) 
_{Mar}
(98) 
_{Apr}
(125) 
_{May}
(54) 
_{Jun}
(119) 
_{Jul}
(60) 
_{Aug}
(111) 
_{Sep}
(118) 
_{Oct}
(125) 
_{Nov}
(119) 
_{Dec}
(94) 
2009 
_{Jan}
(109) 
_{Feb}
(38) 
_{Mar}
(93) 
_{Apr}
(88) 
_{May}
(29) 
_{Jun}
(57) 
_{Jul}
(53) 
_{Aug}
(48) 
_{Sep}
(68) 
_{Oct}
(151) 
_{Nov}
(23) 
_{Dec}
(35) 
2010 
_{Jan}
(84) 
_{Feb}
(60) 
_{Mar}
(184) 
_{Apr}
(112) 
_{May}
(60) 
_{Jun}
(90) 
_{Jul}
(23) 
_{Aug}
(70) 
_{Sep}
(119) 
_{Oct}
(27) 
_{Nov}
(47) 
_{Dec}
(54) 
2011 
_{Jan}
(22) 
_{Feb}
(19) 
_{Mar}
(92) 
_{Apr}
(93) 
_{May}
(35) 
_{Jun}
(91) 
_{Jul}
(32) 
_{Aug}
(61) 
_{Sep}
(7) 
_{Oct}
(69) 
_{Nov}
(81) 
_{Dec}
(23) 
2012 
_{Jan}
(64) 
_{Feb}
(95) 
_{Mar}
(35) 
_{Apr}
(36) 
_{May}
(63) 
_{Jun}
(98) 
_{Jul}
(70) 
_{Aug}
(171) 
_{Sep}
(149) 
_{Oct}
(64) 
_{Nov}
(67) 
_{Dec}
(126) 
2013 
_{Jan}
(108) 
_{Feb}
(104) 
_{Mar}
(171) 
_{Apr}
(133) 
_{May}
(108) 
_{Jun}
(100) 
_{Jul}
(93) 
_{Aug}
(126) 
_{Sep}
(74) 
_{Oct}
(59) 
_{Nov}
(145) 
_{Dec}
(93) 
2014 
_{Jan}
(38) 
_{Feb}
(45) 
_{Mar}
(26) 
_{Apr}
(41) 
_{May}
(125) 
_{Jun}
(70) 
_{Jul}
(61) 
_{Aug}
(66) 
_{Sep}
(60) 
_{Oct}
(110) 
_{Nov}
(27) 
_{Dec}
(30) 
2015 
_{Jan}
(43) 
_{Feb}
(67) 
_{Mar}
(71) 
_{Apr}
(92) 
_{May}
(39) 
_{Jun}
(15) 
_{Jul}
(46) 
_{Aug}

_{Sep}

_{Oct}

_{Nov}

_{Dec}

S  M  T  W  T  F  S 






1
(1) 
2
(5) 
3

4
(2) 
5

6
(16) 
7
(3) 
8
(12) 
9
(2) 
10
(2) 
11

12

13

14
(1) 
15
(7) 
16

17

18
(11) 
19
(6) 
20
(2) 
21
(12) 
22
(4) 
23

24

25
(3) 
26
(3) 
27
(7) 
28
(5) 


From: Roy Stogner <roystgnr@ic...>  20130228 23:36:06

Sorry about the delays, but I'm finally finding some time to work through old emails now. On Thu, 21 Feb 2013, Manav Bhatia wrote: > I am now sure if I understand the concept of "elem_fixed_solution". That's my fault. The FEMSystem time integration APIs were "prematurely optimized", readability really suffered, and we're going to have to change the APIs soon anyway to enable stabilized nonlinear mass matrices. At the moment, in mass_residual(): elem_fixed_solution is equal to some fixed (integratordependent) value of u during your timestep elem_solution is equal to du/dt at the time of evaluation. > [M] = \int_\Omega w^T \du/dt d\Omega + \int_\Omega L(w)^T \tau du/dt d\Omega > > How does the concept of elem_fixed_solution affect this situation with a mass matrix. Assuming your L(w) is really L(u)(w) (e.g. due to advection velocity depending on u) then you'll want to use elem_fixed_solution as the input for that, or you'll want to yell at me to fix the API right now so you can start writing to something saner.  Roy 
From: Roy Stogner <roystgnr@ic...>  20130228 22:54:56

On Thu, 28 Feb 2013, Manav Bhatia wrote: > I wanted to run this by the developers to see if this might be an > unintended behavior: > > Lines 419 in newton_solver.cpp sets the steplength to 1.0, line 420 makes > a call to line_search which returns the value of step length. However, > there is no lvalue speficied to the call to line search, so that steplength > remain 1.0. > > I think the call to line_search should be > > Real steplength = this>line_search(.....); This is embarrassing. As best as I can tell from "svn blame", I committed that bug (a change which ironically enough was being used for debugging purposes) along with the change set I'd been debugging. It looks like we are properly modifying the newton_iterate in line_search itself, though, and the later norm_delta multiplications and divisions by step_length cancel out, so this is only ruining the norm_delta printing in verbose mode? I'll revert it immediately; thank you for catching this!  Roy 
From: Roy Stogner <roystgnr@ic...>  20130228 20:10:13

On Thu, 21 Feb 2013, Manav Bhatia wrote: > What is throwing me off is that the class documentation of Euler2_Solver says: > > * Euler solves u' = f(theta*u_new + (1theta)*u_old), > * Euler2 solves u' = theta*f(u_new) + (1theta)*f(u_old) > > which seems to imply that no mass matrix is accounted for. This is a flaw in the documentation; I'll fix it ASAP. > Also, while going through the code, I came across the following in euler_solver.C (lines 124129) > > // We do a trick here to avoid using a non1 > // elem_solution_derivative: > context.elem_jacobian *= 1.0; > jacobian_computed = _system.mass_residual(jacobian_computed, context) && > jacobian_computed; > context.elem_jacobian *= 1.0; > > > What is the need to multiply the elem_jacobian by 1.0 both before > and after the function call? Wouldn't the call to mass_residual > overwrite the elem_jacobian anyways? The call to mass_residual (or anything else) is supposed to *augment* the Jacobian, rather than writing to a separate jacobian that would then just be added to the original. This was my bright idea for minimizing L1 cache use and is one of the "premature optimization" examples I was talking about. A problem with that is, if mass_residual is computing a jacobian w.r.t. elem_solution, but elem_solution is storing something whose derivative w.r.t. the real solution isn't 1.0, then the user has to remember to multiply those jacobian terms by elem_solution_derivative to take the chain rule into account. The way we do scaling with deltat and 1.0 bits was intended to try and avoid this complication for most codes.  Roy 
From: Manav Bhatia <bhatiamanav@gm...>  20130228 19:26:12

Hi, I wanted to run this by the developers to see if this might be an unintended behavior: Lines 419 in newton_solver.cpp sets the steplength to 1.0, line 420 makes a call to line_search which returns the value of step length. However, there is no lvalue speficied to the call to line search, so that steplength remain 1.0. I think the call to line_search should be Real steplength = this>line_search(.....); Please comment. Thanks, Manav 
From: John Peterson <jwpeterson@gm...>  20130228 00:06:48

On Feb 27, 2013, at 5:24 PM, Carlos Acosta <bdm133@...> wrote: > Thank you very much. That was actually very helpful. I had to save the mesh as quad4 before importing it to libmesh. > > One thing that I found though is that libmesh is not reading the number of local nodes "n_local_nodes()=0". there for the number of degrees of freedom is zero. > I am printing mesh info reading the .e file vs using MeshTools::Generation::build_square (...) and this is what I get. > > *****EXODUS OUTPUT******* > Mesh Information: > mesh_dimension()=2 > spatial_dimension()=3 > n_nodes()=121 > n_local_nodes()=0 > n_elem()=100 > n_local_elem()=0 > n_active_elem()=100 > n_subdomains()=1 > n_partitions()=1 > n_processors()=1 > n_threads()=1 > processor_id()=0 > > EquationSystems > n_systems()=1 > System #0, "Poisson" > Type "LinearImplicit" > Variables="u" > Finite Element Types="LAGRANGE" > Approximation Orders="FIRST" > n_dofs()=0 > n_local_dofs()=0 > n_constrained_dofs()=0 > n_local_constrained_dofs()=0 > n_vectors()=1 > n_matrices()=1 > DofMap Sparsity > Average OnProcessor Bandwidth <= 0 > Average OffProcessor Bandwidth <= 0 > Maximum OnProcessor Bandwidth <= 0 > Maximum OffProcessor Bandwidth <= 0 > DofMap Constraints > Number of DoF Constraints = 0 > ******************* > > vs. > > ****************mesh generation output****** > Mesh Information: > mesh_dimension()=2 > spatial_dimension()=3 > n_nodes()=256 > n_local_nodes()=256 > n_elem()=225 > n_local_elem()=225 > n_active_elem()=225 > n_subdomains()=1 > n_partitions()=1 > n_processors()=1 > n_threads()=1 > processor_id()=0 > > EquationSystems > n_systems()=1 > System #0, "Poisson" > Type "LinearImplicit" > Variables="u" > Finite Element Types="LAGRANGE" > Approximation Orders="FIRST" > n_dofs()=256 > n_local_dofs()=256 > n_constrained_dofs()=0 > n_local_constrained_dofs()=0 > n_vectors()=1 > n_matrices()=1 > DofMap Sparsity > Average OnProcessor Bandwidth <= 8.26562 > Average OffProcessor Bandwidth <= 0 > Maximum OnProcessor Bandwidth <= 9 > Maximum OffProcessor Bandwidth <= 0 > DofMap Constraints > Number of DoF Constraints = 0 > ************************************************ > > Is there an extra step I need to do to read the exodus mesh??? Call mesh.prepare_for_use() maybe? > > Thank you in advance. > > > > > On Wed, Feb 27, 2013 at 1:06 PM, John Peterson <jwpeterson@...> wrote: >> >> >> >> >> On Feb 27, 2013, at 12:23 PM, "Paul T. Bauman" <ptbauman@...> wrote: >> >> > On Wed, Feb 27, 2013 at 11:51 AM, Carlos Acosta <bdm133@...> wrote: >> > >> >> but when I import an exodus mesh I get errors >> >> recognizing the exodus format. >> >> It says: >> >> >> >> >> >> ERROR! Unrecognized element type_str: SHELL4 >> > >> > Hmm, this is strange because Exodus reader should understand SHELL4. >> >> But it doesn't quote work right: shell4's actually have six "sides" and so some of our side id mapping stuff doesn't work. Definitely stick to Quad4's. I think I will patch this to go back to making shell4 throw an error, maybe something like "please export quad4 from cubit instead"... >> >> >> >> > Could >> > you rerun with METHOD=dbg and see if an assert is tripped? How big is your >> > mesh? >> >> It's possible he's using an old version, maybe a sourceforge SVN version? >> >> >> >  >> > Everyone hates slow websites. So do we. >> > Make your web apps faster with AppDynamics >> > Download AppDynamics Lite for free today: >> > http://p.sf.net/sfu/appdyn_d2d_feb >> > _______________________________________________ >> > Libmeshusers mailing list >> > Libmeshusers@... >> > https://lists.sourceforge.net/lists/listinfo/libmeshusers > > > >  > > Carlos Acosta Berlinghieri > Graduate Research Assistant > Department of Biomedical Engineering Computational Center > University of Texas at San Antonio > (210)4584680 (lab phone) 
From: Carlos Acosta <bdm133@my...>  20130227 23:25:24

Thank you very much. That was actually very helpful. I had to save the mesh as quad4 before importing it to libmesh. One thing that I found though is that libmesh is not reading the number of local nodes "n_local_nodes()=0". there for the number of degrees of freedom is zero. I am printing mesh info reading the .e file vs using MeshTools::Generation::build_square (...) and this is what I get. *****EXODUS OUTPUT******* Mesh Information: mesh_dimension()=2 spatial_dimension()=3 n_nodes()=121 n_local_nodes()=0 n_elem()=100 n_local_elem()=0 n_active_elem()=100 n_subdomains()=1 n_partitions()=1 n_processors()=1 n_threads()=1 processor_id()=0 EquationSystems n_systems()=1 System #0, "Poisson" Type "LinearImplicit" Variables="u" Finite Element Types="LAGRANGE" Approximation Orders="FIRST" n_dofs()=0 n_local_dofs()=0 n_constrained_dofs()=0 n_local_constrained_dofs()=0 n_vectors()=1 n_matrices()=1 DofMap Sparsity Average OnProcessor Bandwidth <= 0 Average OffProcessor Bandwidth <= 0 Maximum OnProcessor Bandwidth <= 0 Maximum OffProcessor Bandwidth <= 0 DofMap Constraints Number of DoF Constraints = 0 ******************* vs. ****************mesh generation output****** Mesh Information: mesh_dimension()=2 spatial_dimension()=3 n_nodes()=256 n_local_nodes()=256 n_elem()=225 n_local_elem()=225 n_active_elem()=225 n_subdomains()=1 n_partitions()=1 n_processors()=1 n_threads()=1 processor_id()=0 EquationSystems n_systems()=1 System #0, "Poisson" Type "LinearImplicit" Variables="u" Finite Element Types="LAGRANGE" Approximation Orders="FIRST" n_dofs()=256 n_local_dofs()=256 n_constrained_dofs()=0 n_local_constrained_dofs()=0 n_vectors()=1 n_matrices()=1 DofMap Sparsity Average OnProcessor Bandwidth <= 8.26562 Average OffProcessor Bandwidth <= 0 Maximum OnProcessor Bandwidth <= 9 Maximum OffProcessor Bandwidth <= 0 DofMap Constraints Number of DoF Constraints = 0 ************************************************ Is there an extra step I need to do to read the exodus mesh??? Thank you in advance. On Wed, Feb 27, 2013 at 1:06 PM, John Peterson <jwpeterson@...> wrote: > > > > > On Feb 27, 2013, at 12:23 PM, "Paul T. Bauman" <ptbauman@...> wrote: > > > On Wed, Feb 27, 2013 at 11:51 AM, Carlos Acosta <bdm133@...> > wrote: > > > >> but when I import an exodus mesh I get errors > >> recognizing the exodus format. > >> It says: > >> > >> > >> ERROR! Unrecognized element type_str: SHELL4 > > > > Hmm, this is strange because Exodus reader should understand SHELL4. > > But it doesn't quote work right: shell4's actually have six "sides" and so > some of our side id mapping stuff doesn't work. Definitely stick to > Quad4's. I think I will patch this to go back to making shell4 throw an > error, maybe something like "please export quad4 from cubit instead"... > > > > > Could > > you rerun with METHOD=dbg and see if an assert is tripped? How big is > your > > mesh? > > It's possible he's using an old version, maybe a sourceforge SVN version? > > > > >  > > Everyone hates slow websites. So do we. > > Make your web apps faster with AppDynamics > > Download AppDynamics Lite for free today: > > http://p.sf.net/sfu/appdyn_d2d_feb > > _______________________________________________ > > Libmeshusers mailing list > > Libmeshusers@... > > https://lists.sourceforge.net/lists/listinfo/libmeshusers >  Carlos Acosta Berlinghieri Graduate Research Assistant Department of Biomedical Engineering Computational Center University of Texas at San Antonio (210)4584680 (lab phone) 
From: John Peterson <jwpeterson@gm...>  20130227 19:13:02

On Feb 27, 2013, at 12:23 PM, "Paul T. Bauman" <ptbauman@...> wrote: > On Wed, Feb 27, 2013 at 11:51 AM, Carlos Acosta <bdm133@...> wrote: > >> but when I import an exodus mesh I get errors >> recognizing the exodus format. >> It says: >> >> >> ERROR! Unrecognized element type_str: SHELL4 > > Hmm, this is strange because Exodus reader should understand SHELL4. But it doesn't quote work right: shell4's actually have six "sides" and so some of our side id mapping stuff doesn't work. Definitely stick to Quad4's. I think I will patch this to go back to making shell4 throw an error, maybe something like "please export quad4 from cubit instead"... > Could > you rerun with METHOD=dbg and see if an assert is tripped? How big is your > mesh? It's possible he's using an old version, maybe a sourceforge SVN version? >  > Everyone hates slow websites. So do we. > Make your web apps faster with AppDynamics > Download AppDynamics Lite for free today: > http://p.sf.net/sfu/appdyn_d2d_feb > _______________________________________________ > Libmeshusers mailing list > Libmeshusers@... > https://lists.sourceforge.net/lists/listinfo/libmeshusers 
From: David Knezevic <dknezevic@se...>  20130227 19:05:57

Is there a difference in Cubit between QUAD4 and SHELL4? On 02/27/2013 01:50 PM, John Peterson wrote: > Yes you need to export QUAD4 from cubit instead of SHELL4 which is the default. > >  > John > > On Feb 27, 2013, at 11:51 AM, Carlos Acosta <bdm133@...> wrote: > >> Hello. >> >> I am trying to run one of the first examples. The one solving a Poisson >> problem. >> The example runs great when the mesh is generated from Libmesh >> "MeshTools::Generation::..." but when I import an exodus mesh I get errors >> recognizing the exodus format. >> It says: >> >> >> ERROR! Unrecognized element type_str: SHELL4 >> Stack frames: 6 >> 0: libMesh::print_trace(std::ostream&) >> 1: libMesh::ExodusII_IO_Helper::ElementMaps::assign_conversion(std::string) >> 2: libMesh::ExodusII_IO::read(std::string const&) >> 3: main >> 4: __libc_start_main >> 5: ./exemplo_1.2opt() [0x422049] >> [0] src/mesh/exodusII_io_helper.C, line 1596, compiled Dec 3 2012 at >> 19:00:45 >> terminate called after throwing an instance of 'libMesh::LogicError' >> what(): Error in libMesh internal logic >> [cacostadesktop:16352] *** Process received signal *** >> [cacostadesktop:16352] Signal: Aborted (6) >> [cacostadesktop:16352] Signal code: (6) >> [cacostadesktop:16352] [ 0] /lib/libpthread.so.0(+0xf8f0) [0x7f9b9aaa78f0] >> [cacostadesktop:16352] [ 1] /lib/libc.so.6(gsignal+0x35) [0x7f9b99d99b25] >> [cacostadesktop:16352] [ 2] /lib/libc.so.6(abort+0x180) [0x7f9b99d9d670] >> [cacostadesktop:16352] [ 3] >> /usr/lib/libstdc++.so.6(_ZN9__gnu_cxx27__verbose_terminate_handlerEv+0x115) >> [0x7f9b9d51d8e5] >> [cacostadesktop:16352] [ 4] /usr/lib/libstdc++.so.6(+0xcad16) >> [0x7f9b9d51bd16] >> [cacostadesktop:16352] [ 5] /usr/lib/libstdc++.so.6(+0xcad43) >> [0x7f9b9d51bd43] >> [cacostadesktop:16352] [ 6] /usr/lib/libstdc++.so.6(+0xcae3e) >> [0x7f9b9d51be3e] >> [cacostadesktop:16352] [ 7] >> /usr/share/libmesh/lib/x86_64unknownlinuxgnu_opt/libmesh.so(_ZN7libMesh18ExodusII_IO_Helper11ElementMaps17assign_conversionESs+0x4f0) >> [0x7f9ba2cf1630] >> [cacostadesktop:16352] [ 8] >> /usr/share/libmesh/lib/x86_64unknownlinuxgnu_opt/libmesh.so(_ZN7libMesh11ExodusII_IO4readERKSs+0x31c) >> [0x7f9ba2ceca3c] >> [cacostadesktop:16352] [ 9] ./exemplo_1.2opt(main+0x19f) [0x4269af] >> [cacostadesktop:16352] [10] /lib/libc.so.6(__libc_start_main+0xfd) >> [0x7f9b99d84c4d] >> [cacostadesktop:16352] [11] ./exemplo_1.2opt() [0x422049] >> [cacostadesktop:16352] *** End of error message *** >> Aborted >> >> ***** >> Any clues about what might be missing. >> >> >>  >> >> Carlos Acosta Berlinghieri >> Graduate Research Assistant >> Department of Biomedical Engineering Computational Center >> University of Texas at San Antonio >> (210)4584680 (lab phone) >>  >> Everyone hates slow websites. So do we. >> Make your web apps faster with AppDynamics >> Download AppDynamics Lite for free today: >> http://p.sf.net/sfu/appdyn_d2d_feb >> _______________________________________________ >> Libmeshusers mailing list >> Libmeshusers@... >> https://lists.sourceforge.net/lists/listinfo/libmeshusers >  > Everyone hates slow websites. So do we. > Make your web apps faster with AppDynamics > Download AppDynamics Lite for free today: > http://p.sf.net/sfu/appdyn_d2d_feb > _______________________________________________ > Libmeshusers mailing list > Libmeshusers@... > https://lists.sourceforge.net/lists/listinfo/libmeshusers 
From: Paul T. Bauman <ptbauman@gm...>  20130227 18:53:56

On Wed, Feb 27, 2013 at 11:51 AM, Carlos Acosta <bdm133@...> wrote: > but when I import an exodus mesh I get errors > recognizing the exodus format. > It says: > > > ERROR! Unrecognized element type_str: SHELL4 > Hmm, this is strange because Exodus reader should understand SHELL4. Could you rerun with METHOD=dbg and see if an assert is tripped? How big is your mesh? 
From: John Peterson <jwpeterson@gm...>  20130227 18:51:17

Yes you need to export QUAD4 from cubit instead of SHELL4 which is the default.  John On Feb 27, 2013, at 11:51 AM, Carlos Acosta <bdm133@...> wrote: > Hello. > > I am trying to run one of the first examples. The one solving a Poisson > problem. > The example runs great when the mesh is generated from Libmesh > "MeshTools::Generation::..." but when I import an exodus mesh I get errors > recognizing the exodus format. > It says: > > > ERROR! Unrecognized element type_str: SHELL4 > Stack frames: 6 > 0: libMesh::print_trace(std::ostream&) > 1: libMesh::ExodusII_IO_Helper::ElementMaps::assign_conversion(std::string) > 2: libMesh::ExodusII_IO::read(std::string const&) > 3: main > 4: __libc_start_main > 5: ./exemplo_1.2opt() [0x422049] > [0] src/mesh/exodusII_io_helper.C, line 1596, compiled Dec 3 2012 at > 19:00:45 > terminate called after throwing an instance of 'libMesh::LogicError' > what(): Error in libMesh internal logic > [cacostadesktop:16352] *** Process received signal *** > [cacostadesktop:16352] Signal: Aborted (6) > [cacostadesktop:16352] Signal code: (6) > [cacostadesktop:16352] [ 0] /lib/libpthread.so.0(+0xf8f0) [0x7f9b9aaa78f0] > [cacostadesktop:16352] [ 1] /lib/libc.so.6(gsignal+0x35) [0x7f9b99d99b25] > [cacostadesktop:16352] [ 2] /lib/libc.so.6(abort+0x180) [0x7f9b99d9d670] > [cacostadesktop:16352] [ 3] > /usr/lib/libstdc++.so.6(_ZN9__gnu_cxx27__verbose_terminate_handlerEv+0x115) > [0x7f9b9d51d8e5] > [cacostadesktop:16352] [ 4] /usr/lib/libstdc++.so.6(+0xcad16) > [0x7f9b9d51bd16] > [cacostadesktop:16352] [ 5] /usr/lib/libstdc++.so.6(+0xcad43) > [0x7f9b9d51bd43] > [cacostadesktop:16352] [ 6] /usr/lib/libstdc++.so.6(+0xcae3e) > [0x7f9b9d51be3e] > [cacostadesktop:16352] [ 7] > /usr/share/libmesh/lib/x86_64unknownlinuxgnu_opt/libmesh.so(_ZN7libMesh18ExodusII_IO_Helper11ElementMaps17assign_conversionESs+0x4f0) > [0x7f9ba2cf1630] > [cacostadesktop:16352] [ 8] > /usr/share/libmesh/lib/x86_64unknownlinuxgnu_opt/libmesh.so(_ZN7libMesh11ExodusII_IO4readERKSs+0x31c) > [0x7f9ba2ceca3c] > [cacostadesktop:16352] [ 9] ./exemplo_1.2opt(main+0x19f) [0x4269af] > [cacostadesktop:16352] [10] /lib/libc.so.6(__libc_start_main+0xfd) > [0x7f9b99d84c4d] > [cacostadesktop:16352] [11] ./exemplo_1.2opt() [0x422049] > [cacostadesktop:16352] *** End of error message *** > Aborted > > ***** > Any clues about what might be missing. > > >  > > Carlos Acosta Berlinghieri > Graduate Research Assistant > Department of Biomedical Engineering Computational Center > University of Texas at San Antonio > (210)4584680 (lab phone) >  > Everyone hates slow websites. So do we. > Make your web apps faster with AppDynamics > Download AppDynamics Lite for free today: > http://p.sf.net/sfu/appdyn_d2d_feb > _______________________________________________ > Libmeshusers mailing list > Libmeshusers@... > https://lists.sourceforge.net/lists/listinfo/libmeshusers 
From: Carlos Acosta <bdm133@my...>  20130227 17:51:36

Hello. I am trying to run one of the first examples. The one solving a Poisson problem. The example runs great when the mesh is generated from Libmesh "MeshTools::Generation::..." but when I import an exodus mesh I get errors recognizing the exodus format. It says: ERROR! Unrecognized element type_str: SHELL4 Stack frames: 6 0: libMesh::print_trace(std::ostream&) 1: libMesh::ExodusII_IO_Helper::ElementMaps::assign_conversion(std::string) 2: libMesh::ExodusII_IO::read(std::string const&) 3: main 4: __libc_start_main 5: ./exemplo_1.2opt() [0x422049] [0] src/mesh/exodusII_io_helper.C, line 1596, compiled Dec 3 2012 at 19:00:45 terminate called after throwing an instance of 'libMesh::LogicError' what(): Error in libMesh internal logic [cacostadesktop:16352] *** Process received signal *** [cacostadesktop:16352] Signal: Aborted (6) [cacostadesktop:16352] Signal code: (6) [cacostadesktop:16352] [ 0] /lib/libpthread.so.0(+0xf8f0) [0x7f9b9aaa78f0] [cacostadesktop:16352] [ 1] /lib/libc.so.6(gsignal+0x35) [0x7f9b99d99b25] [cacostadesktop:16352] [ 2] /lib/libc.so.6(abort+0x180) [0x7f9b99d9d670] [cacostadesktop:16352] [ 3] /usr/lib/libstdc++.so.6(_ZN9__gnu_cxx27__verbose_terminate_handlerEv+0x115) [0x7f9b9d51d8e5] [cacostadesktop:16352] [ 4] /usr/lib/libstdc++.so.6(+0xcad16) [0x7f9b9d51bd16] [cacostadesktop:16352] [ 5] /usr/lib/libstdc++.so.6(+0xcad43) [0x7f9b9d51bd43] [cacostadesktop:16352] [ 6] /usr/lib/libstdc++.so.6(+0xcae3e) [0x7f9b9d51be3e] [cacostadesktop:16352] [ 7] /usr/share/libmesh/lib/x86_64unknownlinuxgnu_opt/libmesh.so(_ZN7libMesh18ExodusII_IO_Helper11ElementMaps17assign_conversionESs+0x4f0) [0x7f9ba2cf1630] [cacostadesktop:16352] [ 8] /usr/share/libmesh/lib/x86_64unknownlinuxgnu_opt/libmesh.so(_ZN7libMesh11ExodusII_IO4readERKSs+0x31c) [0x7f9ba2ceca3c] [cacostadesktop:16352] [ 9] ./exemplo_1.2opt(main+0x19f) [0x4269af] [cacostadesktop:16352] [10] /lib/libc.so.6(__libc_start_main+0xfd) [0x7f9b99d84c4d] [cacostadesktop:16352] [11] ./exemplo_1.2opt() [0x422049] [cacostadesktop:16352] *** End of error message *** Aborted ***** Any clues about what might be missing.  Carlos Acosta Berlinghieri Graduate Research Assistant Department of Biomedical Engineering Computational Center University of Texas at San Antonio (210)4584680 (lab phone) 
From: Cody Permann <codypermann@gm...>  20130227 15:53:27

I have a user trying to configure libMesh on a Linux cluster but it keeps running into problems with TBB. He's not alone as I have also had a few other reports of people who still receive error even after adding the diabletbb flag to the configure line. Are we being overlyagressive when snooping for tbb? I haven't really dived into this much yet but just wanted to see if anyone is aware of a quick fix. The compile errros are attached at the bottom of this thread. Thanks, Cody  Forwarded message  From: Capps, Nathan Allen <ncapps@...> Date: Wed, Feb 27, 2013 at 8:30 AM Subject: RE: Compiling problems To: "Permann, Cody J" <cody.permann@...> Cc: "marie.backman@..." <marie.backman@...> I disabled tbb as you suggested. It seem to produce a similar error. Nathan  *From:* Permann, Cody J [cody.permann@...] *Sent:* Wednesday, February 27, 2013 10:17 AM *To:* Capps, Nathan Allen *Subject:* Re: Compiling problems If you take a look at the "build_libmesh_moose.sh" script. It essentially calls libMesh's configure with a couple of switches we need to build moose. You can simply add more flags to that script, in this case "disabletbb". You can always run "./configure help" in libmesh to see what flags and options are available. It's possible you might have to turn off cppthreads too depending on how your system is setup. Cody On Tue, Feb 26, 2013 at 5:29 PM, Capps, Nathan Allen <ncapps@...> wrote: > Thank you for your input and how do I turn off packages which are not > needed and how do I turn off Tbb > > Nathan > > Sent from my iPhone > > On Feb 26, 2013, at 6:28 PM, "Permann, Cody J" <cody.permann@...> > wrote: > > Well  all of your problems are pointing to an issue with TBB (Intel's > Threaded Building Blocks). It looks like libMesh is detecting it but > something isn't quite right with that package. You might want to just > disable (disabletbb). Yeah  you won't be able to run threaded but > that's probably not much of a concern until you are running jobs on several > thousand processors. > > In general, there are a lot of packages that libMesh looks for and > configures if they are on the system. You don't need them all so you can > either try to figure out why it's not working or you can just turn it off. > > BTW  if you want to send output (which is useful) just send it as plain > text :) > > Cody > > > On Tue, Feb 26, 2013 at 4:10 PM, Capps, Nathan Allen <ncapps@...>wrote: > >> Cody, >> >> I apologize for this being a week later. Here is the output when >> compiling lib mesh. I sent you the whole output in hopes to essentially >> help you help me lol. Also I have cc a Brian's new post doc. She is going >> to be using the MOOSE system and is helping me get MOOSE compiled on are >> local cluster. If there is information you need, just let me know. >> >> Nathan >>  >> *From:* Permann, Cody J [cody.permann@...] >> *Sent:* Wednesday, February 20, 2013 3:58 PM >> *To:* Capps, Nathan Allen >> *Subject:* Re: Compiling problems >> >> Are you using our libmesh build script? Is there more to this >> message? It tells me which package is broken but not why. >> >> Cody >> >> >> On Wed, Feb 20, 2013 at 1:35 PM, Capps, Nathan Allen <ncapps@...>wrote: >> >>> Cody, >>> >>> I have been trying to compile the Moose system on my local cluster, >>> and I have been running into a few errors in libmesh listed below. Also I >>> have the petcs3.3p4. Do you have any clues into what the problem might be? >>> >>> collect2: ld returned 1 exit status >>> make[1]: *** [getpot_parseopt] Error 1 >>> make[1]: Leaving directory `/home/ncapps/projects/trunk/libmesh/build' >>> make: *** [installrecursive] Error 1 >>> >>> Nathan >>> >> >> > 
From: John Peterson <jwpeterson@gm...>  20130226 15:55:14

I think using higher order scalar variables is more or less equivalent to using more than one firstorder scalar. It can be a convenient way of grouping together scalar variables which are logically related in some way.  John On Feb 26, 2013, at 8:53 AM, Manav Bhatia <bhatiamanav@...> wrote: > Hi, > > I am trying to understand the SCALAR variable/fe types. > > The FE class of SCALAR seems to return 0 for the shape functions in all > cases. At the same time, the project solution method/class uses > n_components of a SCALAR variable, which in turn is the order of the scalar > variable. > > So, if scalar variable can have higher orders with constant (0.0) shape > function, in which situations can n_components be useful as opposed to > adding multiple scalar variables to a system? > > Any comments would be apprecited. > > Thanks, > Manav >  > Everyone hates slow websites. So do we. > Make your web apps faster with AppDynamics > Download AppDynamics Lite for free today: > http://p.sf.net/sfu/appdyn_d2d_feb > _______________________________________________ > Libmeshusers mailing list > Libmeshusers@... > https://lists.sourceforge.net/lists/listinfo/libmeshusers 
From: Manav Bhatia <bhatiamanav@gm...>  20130226 14:54:03

Hi, I am trying to understand the SCALAR variable/fe types. The FE class of SCALAR seems to return 0 for the shape functions in all cases. At the same time, the project solution method/class uses n_components of a SCALAR variable, which in turn is the order of the scalar variable. So, if scalar variable can have higher orders with constant (0.0) shape function, in which situations can n_components be useful as opposed to adding multiple scalar variables to a system? Any comments would be apprecited. Thanks, Manav 
From: Robert <libmesh@ro...>  20130226 12:00:52

Hello, after looking into it, I think that the MemorySolutionHistory is only usable for timesteps greater than TOLERANCE (which is 1e6 here): > https://github.com/libMesh/libmesh/blob/master/src/solvers/memory_solution_history.C#L108 Is it possible to modify this class to allow smaller timesteps? Regards, Robert 
From: Subramanya Gautam Sadasiva <ssadasiv@pu...>  20130225 21:07:59

Hi, I was looking through the code and I realized that a clear method is called for the nonlinear_solver irrespective of whether it is called through the a NonlinearImplicitSystem or a TransientNonlinearImplicitSystem. This strikes me as a bit wasteful. Is there are a reason for releasing all the memory allocated to the solver, and reconstructing it? I can override it but it seems to be a quite a tedious way to go about it. Is there an easy workaround for this? I am quite new to c++ so it doesn't strike me as a very easy thing to do. Thanks, Subramanya Sadasiva "But memory eventually fades. Turbulences damp out, internal strains yield to plastic flow, concentration inhomogeneities diffuse to uniformity. Systems tend to subside to very simple states,independent of their specific history." Herbert Callen . Thermodynamics and an Introduction to Thermostatics.  Original Message  From: "Subramanya Gautam Sadasiva" <ssadasiv@...> To: "libmeshusers" <libmeshusers@...> Sent: Monday, February 25, 2013 3:39:27 PM Subject: [Libmeshusers] transient nonlinear implicit system Hi, I am trying to use SNESSetPrefix to set a name to the solver to set options for the solver. I do this by deriving from TransientNonlinearImplicitSystem and I do PetscDMNonlinearSolver<Number>* petsc_solver = new PetscDMNonlinearSolver<Number> (*this); // Attach the R & J calculation object petsc_solver>residual_object = this; petsc_solver>jacobian_object = this; // Attach the bounds calculation object petsc_solver>bounds_object = this; // Get a petsc error . this is an alias for a integer PetscErrorCode ierr; petsc_solver>PetscNonlinearSolver<Number>::init(); // The prefix is ch_solve_. I need to figure out how to // put it into the documentation ierr = SNESSetOptionsPrefix(petsc_solver>_snes,"ch_solve_"); in the constructor. This setting lasts for the first call to solve(). However  when I call the solver the second time, all the options are lost. Am I doing something wrong? Do I need to do something else to prevent the solver from being completely reinitialized. Thanks, Subramanya Sadasiva "But memory eventually fades. Turbulences damp out, internal strains yield to plastic flow, concentration inhomogeneities diffuse to uniformity. Systems tend to subside to very simple states,independent of their specific history." Herbert Callen . Thermodynamics and an Introduction to Thermostatics.  Everyone hates slow websites. So do we. Make your web apps faster with AppDynamics Download AppDynamics Lite for free today: http://p.sf.net/sfu/appdyn_d2d_feb _______________________________________________ Libmeshusers mailing list Libmeshusers@... https://lists.sourceforge.net/lists/listinfo/libmeshusers 
From: Subramanya Gautam Sadasiva <ssadasiv@pu...>  20130225 20:40:32

Hi, I am trying to use SNESSetPrefix to set a name to the solver to set options for the solver. I do this by deriving from TransientNonlinearImplicitSystem and I do PetscDMNonlinearSolver<Number>* petsc_solver = new PetscDMNonlinearSolver<Number> (*this); // Attach the R & J calculation object petsc_solver>residual_object = this; petsc_solver>jacobian_object = this; // Attach the bounds calculation object petsc_solver>bounds_object = this; // Get a petsc error . this is an alias for a integer PetscErrorCode ierr; petsc_solver>PetscNonlinearSolver<Number>::init(); // The prefix is ch_solve_. I need to figure out how to // put it into the documentation ierr = SNESSetOptionsPrefix(petsc_solver>_snes,"ch_solve_"); in the constructor. This setting lasts for the first call to solve(). However  when I call the solver the second time, all the options are lost. Am I doing something wrong? Do I need to do something else to prevent the solver from being completely reinitialized. Thanks, Subramanya Sadasiva "But memory eventually fades. Turbulences damp out, internal strains yield to plastic flow, concentration inhomogeneities diffuse to uniformity. Systems tend to subside to very simple states,independent of their specific history." Herbert Callen . Thermodynamics and an Introduction to Thermostatics. 
From: Manav Bhatia <bhatiamanav@gm...>  20130225 04:23:55

Hi, I am attempting to implement my SUPG solver for Euler equations within libMesh via the FEMSystem API. A couple of questions:  I am assuming that when elem_time_derivative is called, elem_solution in the FEMContext stores the u. Similarly, when the mass_residual is called, elem_solution stores the u_dot and elem_fixed_solution stores the u. Is that correct?  How do I get the element u_dot vector when elem_time_derivative is called? I would appreciate any help. Thanks, Manav 
From: Subramanya Gautam Sadasiva <ssadasiv@pu...>  20130222 19:41:44

Sorry about this mail . I see there is a really stupid error in what I am doing. Subramanya Sadasiva "But memory eventually fades. Turbulences damp out, internal strains yield to plastic flow, concentration inhomogeneities diffuse to uniformity. Systems tend to subside to very simple states,independent of their specific history." Herbert Callen . Thermodynamics and an Introduction to Thermostatics.  Original Message  From: "Subramanya Gautam Sadasiva" <ssadasiv@...> To: "libmeshusers" <libmeshusers@...> Sent: Friday, February 22, 2013 2:22:42 PM Subject: [Libmeshusers] Forcing the solver to use a particular solver. Hi, I want to use the petsc_dm solver without passing usepetscdm at the command line. I did this in the constructor of my class that derives from nonlinearimplicitsystem . AutoPtr<PetscNonlinearSolver<Number> > petsc_solver = AutoPtr<PetscNonlinearSolver<Number> > (new PetscNonlinearSolver<Number> (*this)); // Attach the R & J calculation object petsc_solver>residual_object = this; petsc_solver>jacobian_object = this; // Attach the bounds calculation object petsc_solver>bounds_object = this; // Get a petsc error . this is an alias for a integer PetscErrorCode ierr; petsc_solver>PetscNonlinearSolver<Number>::init(); // The prefix is ch_solve_. I need to figure out how to // put it into the documentation ierr = SNESSetOptionsPrefix(petsc_solver>_snes,"ch_solve_"); // Set from options file for the solve ierr = SNESSetFromOptions(petsc_solver>_snes); // Now put petsc solver back into nonlinear solver AutoPtr nonlinear_solver = petsc_solver; This does not seem to work and gives me an error. Is there a better way to do it than to override the nonlinear_solver? Thanks, Subramanya "But memory eventually fades. Turbulences damp out, internal strains yield to plastic flow, concentration inhomogeneities diffuse to uniformity. Systems tend to subside to very simple states,independent of their specific history." Herbert Callen . Thermodynamics and an Introduction to Thermostatics.  Everyone hates slow websites. So do we. Make your web apps faster with AppDynamics Download AppDynamics Lite for free today: http://p.sf.net/sfu/appdyn_d2d_feb _______________________________________________ Libmeshusers mailing list Libmeshusers@... https://lists.sourceforge.net/lists/listinfo/libmeshusers 
From: Subramanya Gautam Sadasiva <ssadasiv@pu...>  20130222 19:22:52

Hi, I want to use the petsc_dm solver without passing usepetscdm at the command line. I did this in the constructor of my class that derives from nonlinearimplicitsystem . AutoPtr<PetscNonlinearSolver<Number> > petsc_solver = AutoPtr<PetscNonlinearSolver<Number> > (new PetscNonlinearSolver<Number> (*this)); // Attach the R & J calculation object petsc_solver>residual_object = this; petsc_solver>jacobian_object = this; // Attach the bounds calculation object petsc_solver>bounds_object = this; // Get a petsc error . this is an alias for a integer PetscErrorCode ierr; petsc_solver>PetscNonlinearSolver<Number>::init(); // The prefix is ch_solve_. I need to figure out how to // put it into the documentation ierr = SNESSetOptionsPrefix(petsc_solver>_snes,"ch_solve_"); // Set from options file for the solve ierr = SNESSetFromOptions(petsc_solver>_snes); // Now put petsc solver back into nonlinear solver AutoPtr nonlinear_solver = petsc_solver; This does not seem to work and gives me an error. Is there a better way to do it than to override the nonlinear_solver? Thanks, Subramanya "But memory eventually fades. Turbulences damp out, internal strains yield to plastic flow, concentration inhomogeneities diffuse to uniformity. Systems tend to subside to very simple states,independent of their specific history." Herbert Callen . Thermodynamics and an Introduction to Thermostatics. 
From: Manav Bhatia <bhatiamanav@gm...>  20130222 04:51:45

On Feb 21, 2013, at 3:25 PM, Roy Stogner <roystgnr@...> wrote: >> >> Additionally, within the existing framework, how does one account for >> *nonlinear* systems that have a mass matrix >> >> [M] du/dt = f(u) > > Right now, we allow you to request an "elem_fixed_solution" that can > be used for either loworder time integration of nonlinear mass > matrices or for nonlinear stabilization in test functions. We're > actually planning a (backwardsincompatible, unfortunately) change in > this API, to make it less unnecessarily confusing and to allow for > higherorder time integration of nonlinear mass matrices and/or for > nonlinear mass matrices with nonlinear stabilization in test > functions. If you'd like, I could probably push that change up higher > on my todo list. Upon further study, it seems like the elem_fixed_solution stores the element solution in situations where the elem_solution may contain the velocity or delta solution. Is this correct? Manav 
From: Manav Bhatia <bhatiamanav@gm...>  20130222 04:44:27

On Feb 21, 2013, at 3:25 PM, Roy Stogner <roystgnr@...> wrote: > >> Additionally, within the existing framework, how does one account for >> *nonlinear* systems that have a mass matrix >> >> [M] du/dt = f(u) > > Right now, we allow you to request an "elem_fixed_solution" that can > be used for either loworder time integration of nonlinear mass > matrices or for nonlinear stabilization in test functions. We're > actually planning a (backwardsincompatible, unfortunately) change in > this API, to make it less unnecessarily confusing and to allow for > higherorder time integration of nonlinear mass matrices and/or for > nonlinear mass matrices with nonlinear stabilization in test > functions. If you'd like, I could probably push that change up higher > on my todo list. > OK, I think I may have started to understand this a little better now (although still a little confused). Please correct me as you see fit. I see that the method mass_residual is used in the classes Euler_Solver and Euler2_Solver. What is throwing me off is that the class documentation of Euler2_Solver says: * Euler solves u' = f(theta*u_new + (1theta)*u_old), * Euler2 solves u' = theta*f(u_new) + (1theta)*f(u_old) which seems to imply that no mass matrix is accounted for. Are the u and f(u) here the solution after discretization, or before discretization? My current guess is that the u' in these equations should actually be replaced with the product [M] u' , like one would encounter postdiscretization. What is your assessment? Also, while going through the code, I came across the following in euler_solver.C (lines 124129) // We do a trick here to avoid using a non1 // elem_solution_derivative: context.elem_jacobian *= 1.0; jacobian_computed = _system.mass_residual(jacobian_computed, context) && jacobian_computed; context.elem_jacobian *= 1.0; What is the need to multiply the elem_jacobian by 1.0 both before and after the function call? Wouldn't the call to mass_residual overwrite the elem_jacobian anyways? Thanks, Manav 
From: Manav Bhatia <bhatiamanav@gm...>  20130221 21:47:42

I had looked at DirichletBoundary class in the examples, but I had not realized that it works for all interpolation basis. What you describe is very impressive, and I will certaintly look deeper into the class. Thanks! Manav On Thu, Feb 21, 2013 at 4:17 PM, Roy Stogner <roystgnr@...>wrote: > > On Thu, 21 Feb 2013, Manav Bhatia wrote: > > I have a preference for Lagrange interpolation for the structural >> (and some thermal) applications that I am working on, primarily due >> to easy application of Dirichlet boundary conditions. >> > > Have you noticed the relativelynew DirichletBoundary class we have? > Hand that a functor, and it'll do the projections onto nonLagrange > spaces for you. Right now there's still work and testing to be done > for transient Dirichlet conditions (that functor may get evaluated at > time t in places where it should be evaluated at t +/ deltat...) but > for steady Dirichlet conditions it's pretty solid. > > > Perhaps hierarchic polynomials would also offer a similar feature >> (please correct if I am wrong) due to the presence of node and edge >> functions whose dof values can be specified, and bubble functions >> that are zero on the boundaries. >> > > Unless your Dirichlet values are piecewiselinear, they'll have > nonzero dofs corresponding to those bubble functions. What our > DirichletBoundary stuff gives you is interpolation at nodes, then > constrained L2 projection for edgeinterior values, then > furtherconstrained L2 projection for faceinterior values. > > > I am curious about your comment "mapped the Lagrange bases on to >> nodes topologically". Are you implying a software fix of telling the >> dof that it is associated with a node, or some sort of a >> mathematical mapping process? >> > > Software. Suppose you have a quad on the unit square. Isoparametric > cubic Lagrange elements would have nodes at (1/3,0) and (2/3,0), each > with a single degree of freedom. My suggestion would be to do what we > do with hierarchics: keep the single node at (1/2,0), give it two > degrees of freedom, then use some canonical ordering based on node > positions to determine which DoF is which. This breaks code which > assumes that solution(node>dof_number(0,0,**0)) == u(Point(node)), but > mathematically it can still be the same basis being used. >  > Roy 
From: Manav Bhatia <bhatiamanav@gm...>  20130221 21:36:48

Hi Roy, In the FEM code that I have written, I have implemented Newmark Method for first and second order nonlinear systems, with and without mass and damping matrices. While converting secondorder to firstorder systems I create a new sparsity patter data structure where the submatrices inherit their patterns from the original system, or get the additional nonzeros that show up from the Identity matrices. I can look at porting this over to libMesh API, and your idea of dof_coupling may serve the necessary purpose here. I can look into this. I am now sure if I understand the concept of "elem_fixed_solution". For the Euler equations L (u) = 0 L = d/dt + A_i d/d x_i that I have implemented, I get the following equations \int_\Omega w^T \du/dt d\Omega  \int_\Omega dw^T/d x_i A_i u d\Omega + \int_\Lambda w^T flux(u) \Lambda + \int_\Omega L(w)^T \tau (d/dt + A_i d/dx_i)u d\Omega = 0 This gives the mass matrix [M] = \int_\Omega w^T \du/dt d\Omega + \int_\Omega L(w)^T \tau du/dt d\Omega and the rhs function f(u) =  ( \int_\Omega dw^T/d x_i A_i u d\Omega + \int_\Lambda w^T flux(u) \Lambda +\int_\Omega L(w)^T \tau A_i du/dx_i d\Omega ) How does the concept of elem_fixed_solution affect this situation with a mass matrix. Now, for steady solutions, I have successfully time marched by assumed the mass matrix to be identity, which brings me back to du/dt = f (u). But, I am curious about timeaccurate simulations, where this assumption may not be accurate. I would appreciate your comments. Thanks, Manav On Thu, Feb 21, 2013 at 3:25 PM, Roy Stogner <roystgnr@...>wrote: > > On Thu, 21 Feb 2013, Manav Bhatia wrote: > > I would appreciate some help in understanding if the existing System >> classes can handle PDEs (linear/nonlinear) with second order time >> derivative. >> > > IIRC We have one existing system (NewmarkSystem?) designed for second > order time derivatives. I have no idea what its limitations are or > whether it would be suitable for you. > > > I am going through the FEMSystem documentation and it seems >> like these are suitable for a system of the form >> >> du/dt = f(u). >> > > Currently, yes. There's nothing about the APIs that would prevent you > from adding a UnsteadySolver subclass that integrates u''=f(u), and > we'd love to get a patch like that, but we don't have any such > subclass in the library yet. > > Worse: if you wanted to handle even more general secondorder time > derivatives, with mixtures of u'' and u' in the same equation, I think > we'd need some modification to the APIs before we could handle that > naturally. > > > Ofcourse, I can change my second order system from >> >> d^2u/dt^2 = f(u) >> >> to >> >> d/dt {u; \tilde{u}} = [0 I; 0 0] {u; \tilde{u}} + {0; f(u)} >> >> but that changes the dimension of my sparse system by a factor of two, >> > > This is probably the easiest way to get a secondorder system up and > running with FEMSystem. > > > and I am not sure if the FEMSystem takes into account the modification in >> the sparsity pattern. >> > > You do the add_variable() for your auxilliary u' variable in your > system initialization, and the sparsity pattern will be modified to > assume that every residual term corresponding to a u degree of freedom > on an element is coupled to every u' degree of freedom on that > element. This will work; you'll just have unnecessary zeros in your > sparse matrix or you'll have to use our DofMap::_dof_coupling API to > avoid them. > > > Additionally, within the existing framework, how does one account for >> *nonlinear* systems that have a mass matrix >> >> [M] du/dt = f(u) >> > > Right now, we allow you to request an "elem_fixed_solution" that can > be used for either loworder time integration of nonlinear mass > matrices or for nonlinear stabilization in test functions. We're > actually planning a (backwardsincompatible, unfortunately) change in > this API, to make it less unnecessarily confusing and to allow for > higherorder time integration of nonlinear mass matrices and/or for > nonlinear mass matrices with nonlinear stabilization in test > functions. If you'd like, I could probably push that change up higher > on my todo list. > > > Extending this question forward, how about PDEs with second order time >> derivative? >> > > Now that's a good question. I've only done secondorder time > integration by breaking it up into firstorder systems. I'd like to > support direct secondorder integrators, but I'd need help designing > the API so as to make it as broadly applicable as possible. >  > Roy > 
From: Roy Stogner <roystgnr@ic...>  20130221 21:17:58

On Thu, 21 Feb 2013, Manav Bhatia wrote: > I have a preference for Lagrange interpolation for the structural > (and some thermal) applications that I am working on, primarily due > to easy application of Dirichlet boundary conditions. Have you noticed the relativelynew DirichletBoundary class we have? Hand that a functor, and it'll do the projections onto nonLagrange spaces for you. Right now there's still work and testing to be done for transient Dirichlet conditions (that functor may get evaluated at time t in places where it should be evaluated at t +/ deltat...) but for steady Dirichlet conditions it's pretty solid. > Perhaps hierarchic polynomials would also offer a similar feature > (please correct if I am wrong) due to the presence of node and edge > functions whose dof values can be specified, and bubble functions > that are zero on the boundaries. Unless your Dirichlet values are piecewiselinear, they'll have nonzero dofs corresponding to those bubble functions. What our DirichletBoundary stuff gives you is interpolation at nodes, then constrained L2 projection for edgeinterior values, then furtherconstrained L2 projection for faceinterior values. > I am curious about your comment "mapped the Lagrange bases on to > nodes topologically". Are you implying a software fix of telling the > dof that it is associated with a node, or some sort of a > mathematical mapping process? Software. Suppose you have a quad on the unit square. Isoparametric cubic Lagrange elements would have nodes at (1/3,0) and (2/3,0), each with a single degree of freedom. My suggestion would be to do what we do with hierarchics: keep the single node at (1/2,0), give it two degrees of freedom, then use some canonical ordering based on node positions to determine which DoF is which. This breaks code which assumes that solution(node>dof_number(0,0,0)) == u(Point(node)), but mathematically it can still be the same basis being used.  Roy 