From: KIRK, B. (JSC-E. (NASA) <ben...@na...> - 2004-02-24 14:34:47
|
The problem seems to be that gcc 3.0 puts the hash map in the std namespace, while all the other 3.x releases put it in __gnu_cxx. Since 3.1 came out so soon after 3.0 I never ran in to this bug. Thanks! in the file include/mesh/mesh_refinement.h, starting at line 328, change it to: #if defined(HAVE_HASH_MAP) typedef std::hash_multimap<unsigned int, Node*> map_type; #elif defined(HAVE_EXT_HASH_MAP) # if __GNUC__ >= 3 # if __GNUC_MINOR__ == 0 typedef std::hash_multimap<unsigned int, Node*> map_type; # else typedef __gnu_cxx::hash_multimap<unsigned int, Node*> map_type; # endif # else DIE A HORRIBLE DEATH # endif #else typedef std::multimap<unsigned int, Node*> map_type; #endif similarly, in src/mesh/mesh_base.C, starting at line 400, change it to: #if defined(HAVE_HASH_MAP) typedef std::hash_multimap<key_type, val_type> map_type; #elif defined(HAVE_EXT_HASH_MAP) # if __GNUC__ >= 3 # if __GNUC_MINOR__ == 0 typedef std::hash_multimap<key_type, val_type> map_type; # else typedef __gnu_cxx::hash_multimap<key_type, val_type> map_type; # endif # else DIE A HORRIBLE DEATH # endif #else typedef std::multimap<key_type, val_type> map_type; #endif Please let me know if this works for you (I don't have access to gcc-3.0 right now). If it does I'll submit the change to CVS so we don't have this problem in the future. Thanks! -Ben -----Original Message----- From: lib...@li... [mailto:lib...@li...] On Behalf Of Rob van Tol Sent: Tuesday, February 24, 2004 7:17 AM To: lib...@li... Subject: [Libmesh-users] (no subject) Dear libmesh users, I just downloaded the code and I am trying to compile it. I am using Redhat 7 and a GCC-3.0 compiler. I get an error message: ISO C++ forbids declaration of `hash_multimap' with no type. Can anyone tell me how to solve this ?? I have included the configure and make messages below. Thank you in advance, Kind regards, ROB configuring libMesh ------------- --------------------------------------------- checking build system type... i686-pc-linux-gnu checking host system type... i686-pc-linux-gnu checking target system type... i686-pc-linux-gnu checking for g++... g++ checking for C++ compiler default output... a.out checking whether the C++ compiler works... yes checking whether we are cross compiling... no checking for suffix of executables... checking for suffix of object files... o checking whether we are using the GNU C++ compiler... yes checking whether g++ accepts -g... yes checking for gcc... gcc checking whether we are using the GNU C compiler... yes checking whether gcc accepts -g... yes checking for gcc option to accept ANSI C... none needed <<< C++ compiler is gcc-3.0 >>> checking how to run the C++ preprocessor... g++ -E checking for egrep... grep -E checking for ANSI C header files... yes checking for sys/types.h... yes checking for sys/stat.h... yes checking for stdlib.h... yes checking for string.h... yes checking for memory.h... yes checking for strings.h... yes checking for inttypes.h... yes checking for stdint.h... yes checking for unistd.h... yes checking for short int... yes checking size of short int... 2 checking for int... yes checking size of int... 4 checking for long int... yes checking size of long int... 4 checking for float... yes checking size of float... 4 checking for double... yes checking size of double... 8 checking for void *... yes checking size of void *... 4 checking getopt.h usability... yes checking getopt.h presence... yes checking for getopt.h... yes checking whether the compiler implements namespaces... yes checking whether the compiler has locale... yes checking whether the compiler has stringstream... yes checking hash_map usability... no checking hash_map presence... no checking for hash_map... no checking ext/hash_map usability... yes checking ext/hash_map presence... yes checking for ext/hash_map... yes checking hash_set usability... no checking hash_set presence... no checking for hash_set... no checking ext/hash_set usability... yes checking ext/hash_set presence... yes checking for ext/hash_set... yes checking zlib.h usability... yes checking zlib.h presence... yes checking for zlib.h... yes <<< Configuring library with AMR support >>> <<< Configuring library with expensive data structures enabled >>> checking rpc/rpc.h usability... yes checking rpc/rpc.h presence... yes checking for rpc/rpc.h... yes <<< Configuring library with XDR support >>> <<< Configuring library with real number support >>> <<< Configuring library with reference counting support >>> checking for ./contrib/netcdf/lib/i686-pc-linux-gnu/libnetcdf.a... no checking for ./contrib/netcdf/include/netcdf.h... no --------------------------------------------- ----- Configuring for optional packages ----- --------------------------------------------- checking for ./contrib/sfcurves/sfcurves.h... yes <<< Configuring library with SFC support >>> checking for ./contrib/gzstream/gzstream.h... yes <<< Configuring library with gzstreams support >>> checking for ./contrib/tecplot/lib/i686-pc-linux-gnu/tecio.a... yes checking for ./contrib/tecplot/include/TECIO.h... yes <<< Configuring library with Tecplot API support >>> checking for ./contrib/laspack/lastypes.h... yes <<< Configuring library with LASPACK version 1.12.3 support >>> checking for /include/petsc.h... no checking for /include/mpi.h... no checking for /lib/libmpich.a... no checking for ./contrib/metis/Lib/metis.h... yes <<< Configuring library with Metis support >>> checking for doxygen... /usr/bin/doxygen checking for dot... no ---------------------------------------------- --- Done configuring for optional packages --- ---------------------------------------------- checking for perl... /usr/bin/perl configure: creating ./config.status config.status: creating Make.common config.status: creating include/base/libmesh_config.h configure: creating ./config.status config.status: creating Make.common config.status: creating doc/Doxyfile config.status: creating include/base/libmesh_config.h config.status: include/base/libmesh_config.h is unchanged --------------------------------------------- --------- Done Configuring libMesh ---------- --------------------------------------------- [vlr@pc02011 libmesh-0.4.1]$ make Compiling C++ (in optimized mode) src/base/dof_map.C... Compiling C++ (in optimized mode) src/base/dof_map_constraints.C... Compiling C++ (in optimized mode) src/base/dof_object.C... Compiling C++ (in optimized mode) src/base/equation_systems.C... Compiling C++ (in optimized mode) src/base/equation_systems_io.C... Compiling C++ (in optimized mode) src/base/frequency_system.C... Compiling C++ (in optimized mode) src/base/libmesh.C... Compiling C++ (in optimized mode) src/base/newmark_system.C... Compiling C++ (in optimized mode) src/base/node.C... Compiling C++ (in optimized mode) src/base/reference_counted_object.C... Compiling C++ (in optimized mode) src/base/reference_counter.C... Compiling C++ (in optimized mode) src/base/steady_system.C... Compiling C++ (in optimized mode) src/base/system_base.C... Compiling C++ (in optimized mode) src/base/system_base_io.C... Compiling C++ (in optimized mode) src/base/system_base_projection.C... Compiling C++ (in optimized mode) src/base/transient_system.C... Compiling C++ (in optimized mode) src/fe/fe_base.C... Compiling C++ (in optimized mode) src/fe/fe_boundary.C... Compiling C++ (in optimized mode) src/fe/fe.C... Compiling C++ (in optimized mode) src/fe/fe_hierarchic.C... Compiling C++ (in optimized mode) src/fe/fe_hierarchic_shape_1D.C... Compiling C++ (in optimized mode) src/fe/fe_hierarchic_shape_2D.C... Compiling C++ (in optimized mode) src/fe/fe_hierarchic_shape_3D.C... Compiling C++ (in optimized mode) src/fe/fe_interface.C... Compiling C++ (in optimized mode) src/fe/fe_interface_inf_fe.C... Compiling C++ (in optimized mode) src/fe/fe_lagrange.C... Compiling C++ (in optimized mode) src/fe/fe_lagrange_shape_1D.C... Compiling C++ (in optimized mode) src/fe/fe_lagrange_shape_2D.C... Compiling C++ (in optimized mode) src/fe/fe_lagrange_shape_3D.C... Compiling C++ (in optimized mode) src/fe/fe_map.C... Compiling C++ (in optimized mode) src/fe/fe_monomial.C... Compiling C++ (in optimized mode) src/fe/fe_monomial_shape_1D.C... Compiling C++ (in optimized mode) src/fe/fe_monomial_shape_2D.C... Compiling C++ (in optimized mode) src/fe/fe_monomial_shape_3D.C... Compiling C++ (in optimized mode) src/fe/inf_fe_base_radial.C... Compiling C++ (in optimized mode) src/fe/inf_fe_boundary.C... Compiling C++ (in optimized mode) src/fe/inf_fe.C... Compiling C++ (in optimized mode) src/fe/inf_fe_jacobi_20_00_eval.C... Compiling C++ (in optimized mode) src/fe/inf_fe_jacobi_30_00_eval.C... Compiling C++ (in optimized mode) src/fe/inf_fe_lagrange_eval.C... Compiling C++ (in optimized mode) src/fe/inf_fe_legendre_eval.C... Compiling C++ (in optimized mode) src/fe/inf_fe_map.C... Compiling C++ (in optimized mode) src/fe/inf_fe_map_eval.C... Compiling C++ (in optimized mode) src/fe/inf_fe_static.C... Compiling C++ (in optimized mode) src/geom/cell.C... Compiling C++ (in optimized mode) src/geom/cell_hex20.C... Compiling C++ (in optimized mode) src/geom/cell_hex27.C... Compiling C++ (in optimized mode) src/geom/cell_hex8.C... Compiling C++ (in optimized mode) src/geom/cell_hex.C... Compiling C++ (in optimized mode) src/geom/cell_inf.C... Compiling C++ (in optimized mode) src/geom/cell_inf_hex16.C... Compiling C++ (in optimized mode) src/geom/cell_inf_hex18.C... Compiling C++ (in optimized mode) src/geom/cell_inf_hex8.C... Compiling C++ (in optimized mode) src/geom/cell_inf_hex.C... Compiling C++ (in optimized mode) src/geom/cell_inf_prism12.C... Compiling C++ (in optimized mode) src/geom/cell_inf_prism6.C... Compiling C++ (in optimized mode) src/geom/cell_inf_prism.C... Compiling C++ (in optimized mode) src/geom/cell_prism15.C... Compiling C++ (in optimized mode) src/geom/cell_prism18.C... Compiling C++ (in optimized mode) src/geom/cell_prism6.C... Compiling C++ (in optimized mode) src/geom/cell_prism.C... Compiling C++ (in optimized mode) src/geom/cell_pyramid5.C... Compiling C++ (in optimized mode) src/geom/cell_pyramid.C... Compiling C++ (in optimized mode) src/geom/cell_tet10.C... Compiling C++ (in optimized mode) src/geom/cell_tet4.C... Compiling C++ (in optimized mode) src/geom/cell_tet.C... Compiling C++ (in optimized mode) src/geom/edge.C... Compiling C++ (in optimized mode) src/geom/edge_edge2.C... Compiling C++ (in optimized mode) src/geom/edge_edge3.C... Compiling C++ (in optimized mode) src/geom/edge_inf_edge2.C... Compiling C++ (in optimized mode) src/geom/elem.C... Compiling C++ (in optimized mode) src/geom/elem_quality.C... Compiling C++ (in optimized mode) src/geom/elem_refinement.C... In file included from src/geom/elem_refinement.C:30: /home/vlr/libmesh-0.4.1/include/mesh/mesh_refinement.h:332: ISO C++ forbids declaration of `hash_multimap' with no type /home/vlr/libmesh-0.4.1/include/mesh/mesh_refinement.h:332: template-id `hash_multimap<unsigned int, Node*>' used as a declarator /home/vlr/libmesh-0.4.1/include/mesh/mesh_refinement.h:332: parse error before `;' token /home/vlr/libmesh-0.4.1/include/mesh/mesh_refinement.h:340: 'map_type' is used as a type, but is not defined as a type. make: *** [src/geom/elem_refinement.i686-pc-linux-gnu.o] Error 1 dr. R. van Tol Computer Simulations and Foundry Processes tel.:+32 (0)9 2645704 - fax:0032 (0)9 2645848 mailto:rob...@wt... <mailto:rob...@wt...> WTCM-CRIF Belgian Center of the Technological Industry Technologiepark 9, B-9052 Gent-Zwijnaarde (Belgium) tel: +32 (0)9/264 5697 http://www.wtcm.be <http://www.wtcm.be/> WTCM-CRIF "Your Gateway To Innovation". Have a look at the future of materials and manufacturing technology http://techniline.wtcm.be <http://techniline.wtcm.be/> |
From: KIRK, B. (JSC-E. (NASA) <ben...@na...> - 2005-10-31 17:41:43
|
Right... For generic quadrature rules we replace \int_\Omega f(x) d\Omega approximately with a sum over {qp} quadrature points. \sum_{qp} f(x_{qp}) w_{qp} J_{qp} where x_{qp} are the locations of the quadrature points in physical space, w_{qp} is the weight of each quadrature point, and J_{qp} is the Jacobian evaluated at the quadrature point (which is not necessarily constant for higher-order mappings from physical to computational space). As John points out, get_JxW returns the product of the weight and the Jacobian for each quadrature point. If you specify a CONSTANT QGauss rule the weight should be 1 and the JxW value should be the Jacobian that you expect. I'd bet you are getting a default 2-point rule that has two points, each weighted by 1/2. As for your second question, I'm not sure I understand, but I think you are asking how to use different orders of mapping for the elements vs. the unknowns? (I'm guessing because of the popular entropy-production you see in 2D subsonic Euler flow over a cylinder when you use a linear map to the reference elements?) In libMesh the elements are mapped with the "natural" Lagrange basis, and the unknowns are what you request. For example, if you solve a problem with CONSTANT MONOMIALS in 2D on Quad9's you will have a solution that is piecewise constant over elements which are mapped quadratically from physical to computational space. -Ben -----Original Message----- From: lib...@li... [mailto:lib...@li...] On Behalf Of John Peterson Sent: Monday, October 31, 2005 10:01 AM To: jc...@ma... Cc: lib...@li... Subject: [Libmesh-users] (no subject) jc...@ma... writes: > Hi, I am trying to implement a Runge-Kutta Discontinuous Galerkin method > to solve the Euler equations in 1D. > > I try to use the FE classes MONOMIAL and XYZ but I get wrong results for > the jacobian evaluation JxW. > > When we transform > > \int_{x0}^{x1} ... dx = \int_{-1}^{+1} ... J dy > > the jacobian is J = dx/dy = ( x1 - x0 ) / 2. Hi, you might be looking at the value of JxW? This is the Jacobian multiplied by the Gauss quadrature weighting value. -J ------------------------------------------------------- This SF.Net email is sponsored by the JBoss Inc. Get Certified Today * Register for a JBoss Training Course Free Certification Exam for All Training Attendees Through End of 2005 Visit http://www.jboss.com/services/certification for more information _______________________________________________ Libmesh-users mailing list Lib...@li... https://lists.sourceforge.net/lists/listinfo/libmesh-users |
From: li p. <li...@ya...> - 2006-07-18 15:54:59
|
I got negative determinant of deformation gradient earlier. This is the problem, which I'm trying to solve for a long time. It seems that the model reaches its limit to deform itself any more. > I'm not sure if I should ask this question here. This is the right list. > Do I need to implement my own error estimator for my > special problem of a complexer geometric domain? > I'm applying Total Largrange formulation to solve a > nonlinear problem. The Kelly error indicator isn't dependent on the domain. If you have a more complex PDE, however, you may need a better indicator. The Kelly indicator is only really rigorous for linear elements solving the Lagrange equation. In particular, if you're using the error indicator as a stopping criterion rather than just to guide refinement, you definitely need a more rigorous error estimator. > My incremental solution can run more than much more > steps without mesh refinement based on > kelly_error_estimator than with it. This sentence hints at the existence of important additional information without actually providing any. If you'd like help you need to fix that. What happens to prevent your solution from running for as many steps as you want? --- Roy Stogner __________________________________________________ Do You Yahoo!? Tired of spam? Yahoo! Mail has the best spam protection around http://mail.yahoo.com |
From: Roy S. <roy...@ic...> - 2006-07-18 16:30:36
|
On Tue, 18 Jul 2006, li pan wrote: > I got negative determinant of deformation gradient > earlier. This is the problem, which I'm trying to > solve for a long time. It seems that the model reaches > its limit to deform itself any more. You're deforming an adaptively refined mesh? That may be a problem unless you're very careful: the libMesh AMR code assumes in a few places that parent elements coincide with the union of their children, and it assumes in a lot of places that hanging nodes haven't been deformed outside of their parent element sides. If your top level mesh is too coarse there won't be a lot of distortion that you can do safely. Even without AMR, however, you're eventually going to get to a point where some elements are too distorted to use. A good mesh smoother will buy you extra time, but I think the only long term solution is to occasionally remesh the distorted domain from scratch and project the solution from the old mesh to the new. --- Roy Stogner |
From: Kurtis P. <Mer...@un...> - 2007-10-01 14:17:41
|
As a business you have been preapproved to receive 48498 USD TODAY! No hassle at all, completely unsecured. There are no hidden costs or fees. Worried that your credit is less than perfect? Not an issue. Give us a ring, now.. 1.877.292.6894 Turn your dream into a reality. 1.877.292.6894 He set the typewriter down, then rocked it up so he could fish out this new surprise. The side of the mower squalled along the side of the cruiser and took off some paint. Vito Frost |
From: Kirk, B. (JSC-EG) <Ben...@na...> - 2008-10-01 00:45:18
|
By default if you declare three systems you will triple the storage, even if the fe types and everything are identical. It would be pretty straightforward to force the matrix storage to be shared by all three systems, but you will still have three times the DofMap storage (and restriction/prolongation overhead in the case of amr). But yes, if you do nothing special you will get triplicate maps which should stay consistent. One of the more hackish ways to handle the situation would be to only declare one system and use a static counter inside your assemble() function to actually assemble (and thus solve) the different systems. Roy or John?? As for the other issue, system.solution is uised as the initial guess to the iterative siolver, so if you copy the contents between systems you should achieve the desired effect. -Ben ----- Original Message ----- From: Adam Arbree <ar...@cs...> To: lib...@li... <lib...@li...> Sent: Tue Sep 30 18:01:40 2008 Subject: [Libmesh-users] (no subject) I have an application where I need to solve three nearly identical systems. For each system, the matrix entries and the rhs will be slightly different, but the mesh and finite element type and degree will be identical. I have three questions: 1) I would like to reuse the DofMap between each of these systems to save memory (and any other memory that is possible, I am not sure what the overhead per system is); 2) if it is not convenient to avoid duplication of the map memory, I would like to construct the system such that at least the three maps are the same (will this happen automatically?); 3) can I specify an initial guess for the solution when using an iterative solver so that I can use the solution of one system as the initial guess for the another. Thanks you very much Adam Arbree ------------------------------------------------------------------------- This SF.Net email is sponsored by the Moblin Your Move Developer's challenge Build the coolest Linux based applications with Moblin SDK & win great prizes Grand prize is a trip for two to an Open Source event anywhere in the world http://moblin-contest.org/redirect.php?banner_id=100&url=/ _______________________________________________ Libmesh-users mailing list Lib...@li... https://lists.sourceforge.net/lists/listinfo/libmesh-users |
From: Adam A. <ar...@cs...> - 2008-10-01 14:20:56
|
Thanks for the help. The actual systems are slightly different so I am willing to pay for triplicate matrix storage. However, I was hoping to save the DofMap storage if possible. That does not seem easy, so I will have to see how close I come to running out of memory. The information on the initial guess will be useful. Adam ________________________________ From: Kirk, Benjamin (JSC-EG) [mailto:Ben...@na...] Sent: Tue 9/30/2008 8:45 PM To: Adam Arbree; lib...@li... Subject: Re: [Libmesh-users] (no subject) By default if you declare three systems you will triple the storage, even if the fe types and everything are identical. It would be pretty straightforward to force the matrix storage to be shared by all three systems, but you will still have three times the DofMap storage (and restriction/prolongation overhead in the case of amr). But yes, if you do nothing special you will get triplicate maps which should stay consistent. One of the more hackish ways to handle the situation would be to only declare one system and use a static counter inside your assemble() function to actually assemble (and thus solve) the different systems. Roy or John?? As for the other issue, system.solution is uised as the initial guess to the iterative siolver, so if you copy the contents between systems you should achieve the desired effect. -Ben ----- Original Message ----- From: Adam Arbree <ar...@cs...> To: lib...@li... <lib...@li...> Sent: Tue Sep 30 18:01:40 2008 Subject: [Libmesh-users] (no subject) I have an application where I need to solve three nearly identical systems. For each system, the matrix entries and the rhs will be slightly different, but the mesh and finite element type and degree will be identical. I have three questions: 1) I would like to reuse the DofMap between each of these systems to save memory (and any other memory that is possible, I am not sure what the overhead per system is); 2) if it is not convenient to avoid duplication of the map memory, I would like to construct the system such that at least the three maps are the same (will this happen automatically?); 3) can I specify an initial guess for the solution when using an iterative solver so that I can use the solution of one system as the initial guess for the another. Thanks you very much Adam Arbree ------------------------------------------------------------------------- This SF.Net email is sponsored by the Moblin Your Move Developer's challenge Build the coolest Linux based applications with Moblin SDK & win great prizes Grand prize is a trip for two to an Open Source event anywhere in the world http://moblin-contest.org/redirect.php?banner_id=100&url=/ _______________________________________________ Libmesh-users mailing list Lib...@li... https://lists.sourceforge.net/lists/listinfo/libmesh-users |
From: John P. <jwp...@gm...> - 2008-10-01 14:52:21
|
> From: Kirk, Benjamin (JSC-EG) [mailto:Ben...@na...] > Sent: Tue 9/30/2008 8:45 PM > To: Adam Arbree; lib...@li... > Subject: Re: [Libmesh-users] (no subject) > One of the more hackish ways to handle the situation would be to only declare > one system and use a static counter inside your assemble() function to > actually assemble (and thus solve) the different systems. Roy or John?? Static may not even be required. I think a single system with 3 different matrix assembly "branches" is what you want. for (qp) for (i) for (j) switch (system_num) // A variable you define to order the 3 systems case 1: K(i,j) += something; break; case 2: K(i,j) += something else; break; end switch end j end i end qp I've done this in the past to assemble different right-hand side vectors depending on whether I was solving the main system or an auxiliary arclength continuation system of equations. It should be straightforward to extend that to the matrix assembly itself as in the pseudocode above. Of course, this litters the code with switch blocks, but I think the DoF savings would be worth it. You will also need to add extra vectors to the system before initialization to store the extra rhs's and solutions, if say there is a 1-way coupling between the different systems. -- John |
From: Roy S. <roy...@ic...> - 2008-10-01 15:41:07
|
On Wed, 1 Oct 2008, John Peterson wrote: > Static may not even be required. I think a single system with 3 > different matrix assembly "branches" is what you want. John is probably exactly right. In the worst case, if one of your systems has a constant matrix, having to reassemble it over and over again will waste a few CPU cycles... but even then it may be worth it to save the memory cost of the entire matrix. > for (qp) > for (i) > for (j) > switch (system_num) // A variable you define to order the 3 systems > case 1: > K(i,j) += something; break; > case 2: > K(i,j) += something else; break; > end switch > end j > end i > end qp I think putting the switch statement outside the for loops is slightly more efficient, if matrix assembly time ever becomes a noticeable issue (as it did for one app of mine). > I've done this in the past to assemble different right-hand side > vectors depending on whether I was solving the main system or an > auxiliary arclength continuation system of equations. It should be > straightforward to extend that to the matrix assembly itself as in the > pseudocode above. It is straightforward, and I have no excuse for not thinking of it and mentioning it myself. From my first shear thinning flow code: void assemble_stokes (EquationSystems& es, const std::string& system_name) { if (find_vorticity) { assemble_vorticity(es, system_name); return; } assemble_stokes() built the residual and jacobian of the Navier-Stokes system to solve; assemble_vorticity built a linear system to project the vorticity into the same smooth space for plotting. --- Roy |
From: Adam A. <ar...@cs...> - 2008-10-01 17:45:09
|
Thank you all for the details. All of these ideas are useful, but I think for my application I will have to build all three matrices simultaneously. I will describe a little more. I am writing an image renderer and the FEM code simulates light transfer in a material using a diffusion process. The three solutions are the three different colors of light the renderer simulates. The matrix assembly requires a lighting simulation performed by another library. This simulation is significantly more expensive than the integration cost and other elements of the assembly. This library produces the all three color elements together and the only practical choice for me is to cache the three color values need for each matrix element until all the solutions are generated. This storage is close to the size of the matrix so I just create three slightly different copies of the matrix and solve them in parallel. Your serial methods would certainly avoid duplication, but I can't really avoid it. However, because I know I must have significant duplication, I want to decrease that footprint as much as possible. In my initial testing with other libraries, I will likely need grids with as many as 1 million elements. My goal for a libmesh implementation is to have as minimal overhead as possible which is why I was hopeful about avoiding DofMap duplication. That said can you give me a sense of how large the DofMap would be for an AMR application, say relative the matrix storage? Also, if you were going to try as quickly as possible, but were willing to do what was necessary to avoid the duplication, what would be the way to start? I also want to position my expertise, I am very comfortable with the programming and such and much less a FEM person, I just really need this work in the next couple of months. Thanks for all the help. Adam > -----Original Message----- > From: Roy Stogner [mailto:roy...@ic...] > Sent: Wednesday, October 01, 2008 11:41 AM > To: John Peterson > Cc: Adam Arbree; lib...@li... > Subject: Re: [Libmesh-users] (no subject) > > > On Wed, 1 Oct 2008, John Peterson wrote: > > > Static may not even be required. I think a single system with 3 > > different matrix assembly "branches" is what you want. > > John is probably exactly right. In the worst case, if one of your > systems has a constant matrix, having to reassemble it over and over > again will waste a few CPU cycles... but even then it may be worth it > to save the memory cost of the entire matrix. > > > for (qp) > > for (i) > > for (j) > > switch (system_num) // A variable you define to order the 3 systems > > case 1: > > K(i,j) += something; break; > > case 2: > > K(i,j) += something else; break; > > end switch > > end j > > end i > > end qp > > I think putting the switch statement outside the for loops is slightly > more efficient, if matrix assembly time ever becomes a noticeable > issue (as it did for one app of mine). > > > I've done this in the past to assemble different right-hand side > > vectors depending on whether I was solving the main system or an > > auxiliary arclength continuation system of equations. It should be > > straightforward to extend that to the matrix assembly itself as in the > > pseudocode above. > > It is straightforward, and I have no excuse for not thinking of it > and mentioning it myself. From my first shear thinning flow code: > > void assemble_stokes (EquationSystems& es, > const std::string& system_name) > { > if (find_vorticity) > { > assemble_vorticity(es, system_name); > return; > } > > assemble_stokes() built the residual and jacobian of the Navier-Stokes > system to solve; assemble_vorticity built a linear system to project > the vorticity into the same smooth space for plotting. > --- > Roy |
From: Roy S. <roy...@ic...> - 2008-10-01 18:12:41
|
On Wed, 1 Oct 2008, Adam Arbree wrote: > That said can you give me a sense of how large the DofMap would be > for an AMR application, say relative the matrix storage? Depends - in serial, IIRC the matrix storage is much larger. In parallel, that changes drastically, since PETSc parallelizes the matrix properly but our SerialMesh class is duplicated on every node. Do you need to use adaptive mesh coarsening? The DofObject costs scale with the local mesh size, and if you only need adaptive mesh refinement I think the bugs there have been worked out of our work-in-progress ParallelMesh. > Also, if you were going to try as quickly as possible, but were > willing to do what was necessary to avoid the duplication, what > would be the way to start? If you're looking for "as quickly as possible", rather than something that could be merged back into the library, it wouldn't be too hard to write a patch that only works for your one problem. Rewrite the DofMap to only write out dof indices if the system number is 0, and rewrite the DofMap and DofObject to pull all requests from system 0 information regardless of what system is really requested. But frankly, I wouldn't worry about it until you've benchmarked things. IIRC I've fit a million elements in a SerialMesh into a few GB of memory before. --- Roy |
From: Adam A. <ar...@cs...> - 2008-10-01 18:30:18
|
That is fantastic; exactly what I needed. Thank you again. Adam > -----Original Message----- > From: Roy Stogner [mailto:roy...@ic...] > Sent: Wednesday, October 01, 2008 2:13 PM > To: Adam Arbree > Cc: lib...@li... > Subject: RE: [Libmesh-users] (no subject) > > > > On Wed, 1 Oct 2008, Adam Arbree wrote: > > > That said can you give me a sense of how large the DofMap would be > > for an AMR application, say relative the matrix storage? > > Depends - in serial, IIRC the matrix storage is much larger. In > parallel, that changes drastically, since PETSc parallelizes the > matrix properly but our SerialMesh class is duplicated on every > node. > > Do you need to use adaptive mesh coarsening? The DofObject costs > scale with the local mesh size, and if you only need adaptive mesh > refinement I think the bugs there have been worked out of our > work-in-progress ParallelMesh. > > > Also, if you were going to try as quickly as possible, but were > > willing to do what was necessary to avoid the duplication, what > > would be the way to start? > > If you're looking for "as quickly as possible", rather than something > that could be merged back into the library, it wouldn't be too > hard to write a patch that only works for your one problem. Rewrite > the DofMap to only write out dof indices if the system number is 0, > and rewrite the DofMap and DofObject to pull all requests from system > 0 information regardless of what system is really requested. > > But frankly, I wouldn't worry about it until you've benchmarked > things. IIRC I've fit a million elements in a SerialMesh into a few > GB of memory before. > --- > Roy |