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}
(28) 
_{Aug}

_{Sep}

_{Oct}

_{Nov}

_{Dec}

S  M  T  W  T  F  S 




1
(3) 
2

3
(4) 
4
(7) 
5

6
(3) 
7
(14) 
8
(1) 
9
(1) 
10
(6) 
11
(1) 
12
(4) 
13
(3) 
14
(1) 
15
(8) 
16
(13) 
17
(21) 
18
(6) 
19
(1) 
20
(2) 
21
(16) 
22
(3) 
23
(12) 
24
(2) 
25

26
(7) 
27
(5) 
28
(12) 
29
(15) 
30

31


From: Jed Brown <jed@59...>  20120829 17:52:56

I'll also comment that aggressive caching is pessimal from a modernhardware and energy perspective. Memory and memory bandwidth take an increasing part of the power and acquisition budget. On modern hardware, for operations that can be vectorized, you should usually plan to recompute everything that takes less than 50 flops per scalar value. If you have a computation that takes more than 50 flops to recompute, you may want to store it, but be aware that reloading it may displace more useful cache lines and if you aren't careful about cache locality (e.g. if the element is visited by a thread running on a different socket later), performance results may be very bad and/or not reproducible. Lei, I suggest being very conservative when selecting what may be worth caching. Also, depending on your application, there may be much larger gains to be had by looking elsewhere. On Wed, Aug 29, 2012 at 12:04 PM, Derek Gaston <friedmud@...> wrote: > With MOOSE we have the full gamut of options. The default is to > recompute everything. We also have the option to cache one FE reinit > in the case where you have a regular grid. We also have the option of > caching EVERY FE reinit for every element so that you can reuse them > until the mesh changes. > > I'll say that that last one is not used very often... it eats up SO > much damn memory... especially in the cases where it would be useful > (like on higher order elements in 3D). So it's kind of a > selfdefeating optimization... > > Derek > > Sent from my iPhone > > On Aug 29, 2012, at 5:56 PM, Roy Stogner <roystgnr@...> wrote: > > > > > On Wed, 29 Aug 2012, Lorenzo Alessio Botti wrote: > > > >> You can always store quantities that you need several times but this > >> might eat up a lot of memory. Accessing to memory might be expensive > >> and sometimes it might be even faster to recompute. > > > > This can absolutely be true  on modern processors some problems end > > up, even in the assembly, bottlenecked by memory bandwidth rather than > > CPU speed. At that point anything that you can compute based on data > > that's already in the processor cache is "free". > > > >> The libMesh approach is to recompute everything. > > > > That's the libMesh examples approach, anyway, but that's intended to > > keep the examples simple more than anything else, I thought. When > > profiling real applications I've found that finding ways to cache > > transcendental function evaluations (e.g. by doing nodal quadrature or > > reinterpolation of such terms) can give a decent savings in runtime. > > > > Of course efficiency discussion is influenced by the fact that all the > > main libMesh developers are profiling based on our own app codes: > > implicit solves with intricate constitutive models. If you're doing > > an explicit calculation then it probably makes sense to cache even > > mappings and shape functions; if you're solving a linear problem then > > you could do much better with Rob Kirby's element matrix > > transformation tricks than with plain quadrature; etc. > >  > > Roy > > > > >  > > 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/ > > _______________________________________________ > > Libmeshusers mailing list > > Libmeshusers@... > > https://lists.sourceforge.net/lists/listinfo/libmeshusers > > >  > 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/ > _______________________________________________ > Libmeshusers mailing list > Libmeshusers@... > https://lists.sourceforge.net/lists/listinfo/libmeshusers > 
From: Derek Gaston <friedmud@gm...>  20120829 17:04:57

With MOOSE we have the full gamut of options. The default is to recompute everything. We also have the option to cache one FE reinit in the case where you have a regular grid. We also have the option of caching EVERY FE reinit for every element so that you can reuse them until the mesh changes. I'll say that that last one is not used very often... it eats up SO much damn memory... especially in the cases where it would be useful (like on higher order elements in 3D). So it's kind of a selfdefeating optimization... Derek Sent from my iPhone On Aug 29, 2012, at 5:56 PM, Roy Stogner <roystgnr@...> wrote: > > On Wed, 29 Aug 2012, Lorenzo Alessio Botti wrote: > >> You can always store quantities that you need several times but this >> might eat up a lot of memory. Accessing to memory might be expensive >> and sometimes it might be even faster to recompute. > > This can absolutely be true  on modern processors some problems end > up, even in the assembly, bottlenecked by memory bandwidth rather than > CPU speed. At that point anything that you can compute based on data > that's already in the processor cache is "free". > >> The libMesh approach is to recompute everything. > > That's the libMesh examples approach, anyway, but that's intended to > keep the examples simple more than anything else, I thought. When > profiling real applications I've found that finding ways to cache > transcendental function evaluations (e.g. by doing nodal quadrature or > reinterpolation of such terms) can give a decent savings in runtime. > > Of course efficiency discussion is influenced by the fact that all the > main libMesh developers are profiling based on our own app codes: > implicit solves with intricate constitutive models. If you're doing > an explicit calculation then it probably makes sense to cache even > mappings and shape functions; if you're solving a linear problem then > you could do much better with Rob Kirby's element matrix > transformation tricks than with plain quadrature; etc. >  > Roy > >  > 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/ > _______________________________________________ > Libmeshusers mailing list > Libmeshusers@... > https://lists.sourceforge.net/lists/listinfo/libmeshusers 
From: Roy Stogner <roystgnr@ic...>  20120829 16:56:31

On Wed, 29 Aug 2012, Lorenzo Alessio Botti wrote: > You can always store quantities that you need several times but this > might eat up a lot of memory. Accessing to memory might be expensive > and sometimes it might be even faster to recompute. This can absolutely be true  on modern processors some problems end up, even in the assembly, bottlenecked by memory bandwidth rather than CPU speed. At that point anything that you can compute based on data that's already in the processor cache is "free". > The libMesh approach is to recompute everything. That's the libMesh examples approach, anyway, but that's intended to keep the examples simple more than anything else, I thought. When profiling real applications I've found that finding ways to cache transcendental function evaluations (e.g. by doing nodal quadrature or reinterpolation of such terms) can give a decent savings in runtime. Of course efficiency discussion is influenced by the fact that all the main libMesh developers are profiling based on our own app codes: implicit solves with intricate constitutive models. If you're doing an explicit calculation then it probably makes sense to cache even mappings and shape functions; if you're solving a linear problem then you could do much better with Rob Kirby's element matrix transformation tricks than with plain quadrature; etc.  Roy 
From: Lorenzo Alessio Botti <bottilorenzo@gm...>  20120829 08:48:41

Hi Lei, I'll try to answer to some of your questions since I wrote this dG example. In general it is not easy to identify which is the most expensive portion of a code and libMesh developers make extensive use of profiling tools. Something might be inefficient in your specific application but efficiency concerns imply that you have something to compare with. > > In this example, the Laplace equation is solved using DG with interior penalty method. So the main job is to compute the mass matrix from the element volume integration and surface jump integration. It is done element by element with the help of element iterator. I know all of the elements organized as a quadforest and only the active items need to be looped. I guess you mean stiffness matrix. > > so does this mean when we need to reach to the next active element, const_element_iterator++ need to go through the quadforest and pass several nonactive parent elements then get to the correct element pointer to the next active element? If the answer is true. Is it too expensive to find an element in a solver? > > The other thing is the metrics computing such as Jacobian, xi_x xi_y eta_x eta_y and face normal etc. Do those metrics computed several times using FE::FEMap? If the mesh is fixed in space, we only need to compute them once and store them somewhere. > > Same thing happens when computing quadrature points in physical domain and value of shape function/ its derivation at those points in computational domain. I think those are computed in >  > fe>reinit(elem) >  > Those values do not change during the calculation and should be computed only once. I think they are in the deepest loop and will eat a lot of cpu time. > You can always store quantities that you need several times but this might eat up a lot of memory. Accessing to memory might be expensive and sometimes it might be even faster to recompute. The libMesh approach is to recompute everything. In your application you might store something but, when the problem you need to solve reaches a critical size, you want to store quantities that are really expensive to compute. I think that one important point is time integration. In explicit solvers assembly times dominates and it might be useful to store FE related quantities. In implicit solvers usually the memory is used to store the global system matrix, the preconditioner ecc... > Last one, to compute the flux jump at flux points on the common surface of two cells with different h levels. The coordinates and normals of the flux quadrature points on the both side of element are computed several times if we use a residual based scheme. It does not change at each refinement level. Is there any efficient way to handle those kind of situation: flux jump calculation at a surface with hanging node? When you do hrefinement the calculation of interface quantities is actually more expensive since you need to find matching quadrature points on the sides of neighboring elements. If you use reference frame basis functions starting from a physical frame quadrature point you need to find the corresponding reference frame quadrature point on the neighbor. This requires to invert the reference to physical frame map by means of the Newton method. > > Or those are just for demo only. Those efficiency problem get even worse for a scheme need residual jacobian at each time step. Do we need to use FEMap explicitly and store the metrics for each element/face? So we don't need to compute them more than once. This is just a demo that make sense for the Laplace equation and it shows how to compute the quantities you need. If you want to store quantities you will need to care care of that. > > Thanks a lot for your time. Greetings Lorenzo > > > > Sincerely Yours, > > Lei Shi >  >  > 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/_______________________________________________ > Libmeshdevel mailing list > Libmeshdevel@... > https://lists.sourceforge.net/lists/listinfo/libmeshdevel 
From: 蔡园武 <yuanwucai@gm...>  20120829 06:03:53

Thanks for your detailed explanation. I tried your method and my code works right now! 2012/8/29 Roy Stogner <roystgnr@...>: > > On Wed, 29 Aug 2012, 蔡园武 wrote: > >> I found the definition of Parallel::barrier: >> void libMesh::Parallel::barrier ( const Communicator & >> comm = >> Communicator_World ) [inline] >> and the explanation: >> Pause execution until all processors reach a certain point. > > >> But I don't know what should be 'a certain point'? > > > The barrier() call. No processor gets to continue past the barrier > until all other processors have reached the barrier. > > >> What can I pass as variable 'Communicator_World'? > > > That's a C++ argumentwithadefaultvalue. You don't pass anything > to it unless you want to use a nondefault value. Just use > > Parallel::barrier(); > > and you're done. >  > Roy  Cai Yuanwu 蔡园武 Dept. of Engineering Mechanics, Dalian University of Technology, Dalian 116024, China 
From: Roy Stogner <roystgnr@ic...>  20120829 04:56:19

On Wed, 29 Aug 2012, 蔡园武 wrote: > I found the definition of Parallel::barrier: > void libMesh::Parallel::barrier ( const Communicator & comm = > Communicator_World ) [inline] > and the explanation: > Pause execution until all processors reach a certain point. > But I don't know what should be 'a certain point'? The barrier() call. No processor gets to continue past the barrier until all other processors have reached the barrier. > What can I pass as variable 'Communicator_World'? That's a C++ argumentwithadefaultvalue. You don't pass anything to it unless you want to use a nondefault value. Just use Parallel::barrier(); and you're done.  Roy 
From: Vikram Garg <vikram.garg@gm...>  20120829 02:46:40

Yuanwu, There is a variable of system that you can set to use the same Ke, its called assemble_before_solve, so simply set this variable to false using the equation systems object you have. That will take care of the multiple assembles. To avoid the multiple solves, please use the following steps: 1) // Get the system's linear solver object to set the preconditioner reuse flag LinearSolver<Number> *linear_solver = system.get_linear_solver(); 2) // Set the reuse_preconditioner variable to true linear_solver>reuse_preconditioner(true); Give that a try. On Tue, Aug 28, 2012 at 9:36 PM, 蔡园武 <yuanwucai@...> wrote: > Good. Anything else I should set in es or system in my code to use the > same matrix K for different Fe, besides 'pc_type lu' in command? > Temporarily I assemble K for each Fe. But if I use LU, there's no need > to assemble K for multiple times. > > 2012/8/29 Vikram Garg <vikram.v.garg@...>: > > Well, if you use lu, then there should be no 'residual' provided the > matrix > > is invertible. > > > >  > Cai Yuanwu 蔡园武 > Dept. of Engineering Mechanics, > Dalian University of Technology, > Dalian 116024, China >  Vikram Garg http://users.ices.utexas.edu/~vikram/ http://www.runforindia.org/runners/vikramg 
From: 蔡园武 <yuanwucai@gm...>  20120829 02:44:47

2012/8/27 Roy Stogner <roystgnr@...>: > > On Mon, 27 Aug 2012, 蔡园武 wrote: > >> It seems that the processor_id 0 is writing equation_systems to file >> "Iter_0.gmv", and at the same time, processor_id 1 is trying to open >> the file "Iter_0.gmv". Is it correct? > > > That would be my first guess. > > >> I don't know much about MPI. Could someone explain this error, > > > For serialized file formats, we gather data to processor_id 0 and then > other processors are allowed to continue while 0 writes. > > >> and how to solve it? > > > Try a Parallel::barrier() call after the GMVIO? >  > Roy I found the definition of Parallel::barrier: void libMesh::Parallel::barrier ( const Communicator & comm = Communicator_World ) [inline] and the explanation: Pause execution until all processors reach a certain point. But I don't know what should be 'a certain point'? What can I pass as variable 'Communicator_World'? Thanks!  Cai Yuanwu 蔡园武 Dept. of Engineering Mechanics, Dalian University of Technology, Dalian 116024, China 
From: 蔡园武 <yuanwucai@gm...>  20120829 02:36:18

Good. Anything else I should set in es or system in my code to use the same matrix K for different Fe, besides 'pc_type lu' in command? Temporarily I assemble K for each Fe. But if I use LU, there's no need to assemble K for multiple times. 2012/8/29 Vikram Garg <vikram.v.garg@...>: > Well, if you use lu, then there should be no 'residual' provided the matrix > is invertible. >  Cai Yuanwu 蔡园武 Dept. of Engineering Mechanics, Dalian University of Technology, Dalian 116024, China 
From: Vikram Garg <vikram.garg@gm...>  20120829 02:26:08

Well, if you use lu, then there should be no 'residual' provided the matrix is invertible. On Tue, Aug 28, 2012 at 9:24 PM, 蔡园武 <yuanwucai@...> wrote: > Good idea, man, but I don't think changing preconditioner can control > the residual error. > > 2012/8/29 Vikram Garg <simulationist@...>: > > Yuanwu, > > To set the lu preconditioner type just add the following to > > the command you use to execute: pc_type lu . If you want to use some > other > > kind of preconditioner you can just change lu to whatever type you want > > (ilu, jacobi etc) I believe the default is jacobi. > > > > Warning: If your problem sizes are really large then lu might take a > while. > > Just try it out first and then you can try ilu if its too slow. > > > > Thanks. > > > > > > On Tue, Aug 28, 2012 at 9:02 PM, 蔡园武 <yuanwucai@...> wrote: > >> > >> 2012/8/29 Vikram Garg <simulationist@...>: > >> > Hey Yuanwu, > >> > Can you tell us what kind of linear solver > options > >> > you > >> > are using (the method, preconditioner options passed at runtime) ? > >> > >> I used the default option. (I don't know what's the default?) > >> #EquationSystems > >> n_systems()=1 > >> System #0, "homogenization_elastic_3D_plate" > >> Type "LinearImplicit" > >> Variables="u" "v" "w" > >> Finite Element Types="LAGRANGE" "LAGRANGE" "LAGRANGE" > >> Approximation Orders="FIRST" "FIRST" "FIRST" > >> n_dofs()=648 > >> n_local_dofs()=648 > >> n_constrained_dofs()=198 > >> n_local_constrained_dofs()=198 > >> n_vectors()=13 > >> n_matrices()=1 > >> DofMap Sparsity > >> Average OnProcessor Bandwidth <= 84 > >> Average OffProcessor Bandwidth <= 0 > >> Maximum OnProcessor Bandwidth <= 144 > >> Maximum OffProcessor Bandwidth <= 0 > >> DofMap Constraints > >> Number of DoF Constraints = 198 > >> Average DoF Constraint Length= 1 > >> > >> > If you > >> > are using the same stiffness matrix (K) and just changing the rhs (F), > >> > >> That's just what I'm doing. > >> > >> > then > >> > you might actually be better of just using the LU preconditioner for > one > >> > load and reusing it for the others. > >> > >> How to set it in my code? Can you give me a simple example? > >> Thanks for your answer! > >> > >> > > >> > Thanks. > >> > > >> > On Tue, Aug 28, 2012 at 8:45 PM, 蔡园武 <yuanwucai@...> wrote: > >> >> > >> >> Hi, guys, > >> >> I have a LinearImplicitSystem, solved with a sequential of different > >> >> force vectors (rhs). > >> >> Actually I defined different 'Fe' assemble function using a > 'loadcase' > >> >> indicator. In main function, I set the 'loadcase' value, call > >> >> system.solve(), then es.reinit(), set a new 'loadcase' value, and > >> >> solve() again. > >> >> But I found that for some 'loadcase', the solver converged badly, > like: > >> >> > >> >> loadcase1: Linear solver converged at step: 10736, final residual: > >> >> 2.00331e21 > >> >> loadcase2: Linear solver converged at step: 8549, final residual: > >> >> 2.10685e21 > >> >> loadcase3: Linear solver converged at step: 8, final residual: > >> >> 1.38269e07 > >> >> loadcase4: Linear solver converged at step: 0, final residual: > 0.463112 > >> >> loadcase5: Linear solver converged at step: 0, final residual: > 0.463112 > >> >> loadcase6: Linear solver converged at step: 0, final residual: > 0.463112 > >> >> > >> >> In loadcase3, Linear solver converged at step 8, final resudual is > not > >> >> small enough. It's very strange that the loadcases after this didn't > >> >> run? (converged at step 0?) The results are wrong and unbelievable. > >> >> But if I solve loadcase3 after all the other loadcases, then loadcase > >> >> 4,5,6 will be all right. Don't know why? > >> >> > >> >> Can I mannually control the final resudual tolerance? I did set a > >> >> parameter in es: > >> >> es.parameters.set<unsigned int> ("linear solver maximum > >> >> iterations") = 20000; > >> >> es.parameters.set<Real> ("linear solver tolerance") = TOLERANCE; > >> >> > >> >> Thanks for your help！ > >> >>  > >> >> Cai Yuanwu 蔡园武 > >> >> Dept. of Engineering Mechanics, > >> >> Dalian University of Technology, > >> >> Dalian 116024, China > >> >> > >> >> > >> >> > >> >> >  > >> >> 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/ > >> >> _______________________________________________ > >> >> Libmeshusers mailing list > >> >> Libmeshusers@... > >> >> https://lists.sourceforge.net/lists/listinfo/libmeshusers > >> > > >> > > >> > > >> > > >> >  > >> > Vikram Garg > >> > PhD Candidate > >> > Institute for Computational and Engineering Sciences > >> > The University of Texas at Austin > >> > > >> > http://users.ices.utexas.edu/~vikram/ > >> > http://www.runforindia.org/runners/vikramg > >> > > >> > >> > >> > >>  > >> Cai Yuanwu 蔡园武 > >> Dept. of Engineering Mechanics, > >> Dalian University of Technology, > >> Dalian 116024, China > > > > > > > > > >  > > Vikram Garg > > PhD Candidate > > Institute for Computational and Engineering Sciences > > The University of Texas at Austin > > > > http://users.ices.utexas.edu/~vikram/ > > http://www.runforindia.org/runners/vikramg > > > > > >  > Cai Yuanwu 蔡园武 > Dept. of Engineering Mechanics, > Dalian University of Technology, > Dalian 116024, China >  Vikram Garg http://users.ices.utexas.edu/~vikram/ http://www.runforindia.org/runners/vikramg 
From: 蔡园武 <yuanwucai@gm...>  20120829 02:24:38

Good idea, man, but I don't think changing preconditioner can control the residual error. 2012/8/29 Vikram Garg <simulationist@...>: > Yuanwu, > To set the lu preconditioner type just add the following to > the command you use to execute: pc_type lu . If you want to use some other > kind of preconditioner you can just change lu to whatever type you want > (ilu, jacobi etc) I believe the default is jacobi. > > Warning: If your problem sizes are really large then lu might take a while. > Just try it out first and then you can try ilu if its too slow. > > Thanks. > > > On Tue, Aug 28, 2012 at 9:02 PM, 蔡园武 <yuanwucai@...> wrote: >> >> 2012/8/29 Vikram Garg <simulationist@...>: >> > Hey Yuanwu, >> > Can you tell us what kind of linear solver options >> > you >> > are using (the method, preconditioner options passed at runtime) ? >> >> I used the default option. (I don't know what's the default?) >> #EquationSystems >> n_systems()=1 >> System #0, "homogenization_elastic_3D_plate" >> Type "LinearImplicit" >> Variables="u" "v" "w" >> Finite Element Types="LAGRANGE" "LAGRANGE" "LAGRANGE" >> Approximation Orders="FIRST" "FIRST" "FIRST" >> n_dofs()=648 >> n_local_dofs()=648 >> n_constrained_dofs()=198 >> n_local_constrained_dofs()=198 >> n_vectors()=13 >> n_matrices()=1 >> DofMap Sparsity >> Average OnProcessor Bandwidth <= 84 >> Average OffProcessor Bandwidth <= 0 >> Maximum OnProcessor Bandwidth <= 144 >> Maximum OffProcessor Bandwidth <= 0 >> DofMap Constraints >> Number of DoF Constraints = 198 >> Average DoF Constraint Length= 1 >> >> > If you >> > are using the same stiffness matrix (K) and just changing the rhs (F), >> >> That's just what I'm doing. >> >> > then >> > you might actually be better of just using the LU preconditioner for one >> > load and reusing it for the others. >> >> How to set it in my code? Can you give me a simple example? >> Thanks for your answer! >> >> > >> > Thanks. >> > >> > On Tue, Aug 28, 2012 at 8:45 PM, 蔡园武 <yuanwucai@...> wrote: >> >> >> >> Hi, guys, >> >> I have a LinearImplicitSystem, solved with a sequential of different >> >> force vectors (rhs). >> >> Actually I defined different 'Fe' assemble function using a 'loadcase' >> >> indicator. In main function, I set the 'loadcase' value, call >> >> system.solve(), then es.reinit(), set a new 'loadcase' value, and >> >> solve() again. >> >> But I found that for some 'loadcase', the solver converged badly, like: >> >> >> >> loadcase1: Linear solver converged at step: 10736, final residual: >> >> 2.00331e21 >> >> loadcase2: Linear solver converged at step: 8549, final residual: >> >> 2.10685e21 >> >> loadcase3: Linear solver converged at step: 8, final residual: >> >> 1.38269e07 >> >> loadcase4: Linear solver converged at step: 0, final residual: 0.463112 >> >> loadcase5: Linear solver converged at step: 0, final residual: 0.463112 >> >> loadcase6: Linear solver converged at step: 0, final residual: 0.463112 >> >> >> >> In loadcase3, Linear solver converged at step 8, final resudual is not >> >> small enough. It's very strange that the loadcases after this didn't >> >> run? (converged at step 0?) The results are wrong and unbelievable. >> >> But if I solve loadcase3 after all the other loadcases, then loadcase >> >> 4,5,6 will be all right. Don't know why? >> >> >> >> Can I mannually control the final resudual tolerance? I did set a >> >> parameter in es: >> >> es.parameters.set<unsigned int> ("linear solver maximum >> >> iterations") = 20000; >> >> es.parameters.set<Real> ("linear solver tolerance") = TOLERANCE; >> >> >> >> Thanks for your help！ >> >>  >> >> Cai Yuanwu 蔡园武 >> >> Dept. of Engineering Mechanics, >> >> Dalian University of Technology, >> >> Dalian 116024, China >> >> >> >> >> >> >> >>  >> >> 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/ >> >> _______________________________________________ >> >> Libmeshusers mailing list >> >> Libmeshusers@... >> >> https://lists.sourceforge.net/lists/listinfo/libmeshusers >> > >> > >> > >> > >> >  >> > Vikram Garg >> > PhD Candidate >> > Institute for Computational and Engineering Sciences >> > The University of Texas at Austin >> > >> > http://users.ices.utexas.edu/~vikram/ >> > http://www.runforindia.org/runners/vikramg >> > >> >> >> >>  >> Cai Yuanwu 蔡园武 >> Dept. of Engineering Mechanics, >> Dalian University of Technology, >> Dalian 116024, China > > > > >  > Vikram Garg > PhD Candidate > Institute for Computational and Engineering Sciences > The University of Texas at Austin > > http://users.ices.utexas.edu/~vikram/ > http://www.runforindia.org/runners/vikramg >  Cai Yuanwu 蔡园武 Dept. of Engineering Mechanics, Dalian University of Technology, Dalian 116024, China 
From: Vikram Garg <simulationist@gm...>  20120829 02:07:14

Yuanwu, To set the lu preconditioner type just add the following to the command you use to execute: pc_type lu . If you want to use some other kind of preconditioner you can just change lu to whatever type you want (ilu, jacobi etc) I believe the default is jacobi. Warning: If your problem sizes are really large then lu might take a while. Just try it out first and then you can try ilu if its too slow. Thanks. On Tue, Aug 28, 2012 at 9:02 PM, 蔡园武 <yuanwucai@...> wrote: > 2012/8/29 Vikram Garg <simulationist@...>: > > Hey Yuanwu, > > Can you tell us what kind of linear solver options > you > > are using (the method, preconditioner options passed at runtime) ? > > I used the default option. (I don't know what's the default?) > #EquationSystems > n_systems()=1 > System #0, "homogenization_elastic_3D_plate" > Type "LinearImplicit" > Variables="u" "v" "w" > Finite Element Types="LAGRANGE" "LAGRANGE" "LAGRANGE" > Approximation Orders="FIRST" "FIRST" "FIRST" > n_dofs()=648 > n_local_dofs()=648 > n_constrained_dofs()=198 > n_local_constrained_dofs()=198 > n_vectors()=13 > n_matrices()=1 > DofMap Sparsity > Average OnProcessor Bandwidth <= 84 > Average OffProcessor Bandwidth <= 0 > Maximum OnProcessor Bandwidth <= 144 > Maximum OffProcessor Bandwidth <= 0 > DofMap Constraints > Number of DoF Constraints = 198 > Average DoF Constraint Length= 1 > > > If you > > are using the same stiffness matrix (K) and just changing the rhs (F), > > That's just what I'm doing. > > > then > > you might actually be better of just using the LU preconditioner for one > > load and reusing it for the others. > > How to set it in my code? Can you give me a simple example? > Thanks for your answer! > > > > > Thanks. > > > > On Tue, Aug 28, 2012 at 8:45 PM, 蔡园武 <yuanwucai@...> wrote: > >> > >> Hi, guys, > >> I have a LinearImplicitSystem, solved with a sequential of different > >> force vectors (rhs). > >> Actually I defined different 'Fe' assemble function using a 'loadcase' > >> indicator. In main function, I set the 'loadcase' value, call > >> system.solve(), then es.reinit(), set a new 'loadcase' value, and > >> solve() again. > >> But I found that for some 'loadcase', the solver converged badly, like: > >> > >> loadcase1: Linear solver converged at step: 10736, final residual: > >> 2.00331e21 > >> loadcase2: Linear solver converged at step: 8549, final residual: > >> 2.10685e21 > >> loadcase3: Linear solver converged at step: 8, final residual: > 1.38269e07 > >> loadcase4: Linear solver converged at step: 0, final residual: 0.463112 > >> loadcase5: Linear solver converged at step: 0, final residual: 0.463112 > >> loadcase6: Linear solver converged at step: 0, final residual: 0.463112 > >> > >> In loadcase3, Linear solver converged at step 8, final resudual is not > >> small enough. It's very strange that the loadcases after this didn't > >> run? (converged at step 0?) The results are wrong and unbelievable. > >> But if I solve loadcase3 after all the other loadcases, then loadcase > >> 4,5,6 will be all right. Don't know why? > >> > >> Can I mannually control the final resudual tolerance? I did set a > >> parameter in es: > >> es.parameters.set<unsigned int> ("linear solver maximum > >> iterations") = 20000; > >> es.parameters.set<Real> ("linear solver tolerance") = TOLERANCE; > >> > >> Thanks for your help！ > >>  > >> Cai Yuanwu 蔡园武 > >> Dept. of Engineering Mechanics, > >> Dalian University of Technology, > >> Dalian 116024, China > >> > >> > >> >  > >> 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/ > >> _______________________________________________ > >> Libmeshusers mailing list > >> Libmeshusers@... > >> https://lists.sourceforge.net/lists/listinfo/libmeshusers > > > > > > > > > >  > > Vikram Garg > > PhD Candidate > > Institute for Computational and Engineering Sciences > > The University of Texas at Austin > > > > http://users.ices.utexas.edu/~vikram/ > > http://www.runforindia.org/runners/vikramg > > > > > >  > Cai Yuanwu 蔡园武 > Dept. of Engineering Mechanics, > Dalian University of Technology, > Dalian 116024, China >  Vikram Garg PhD Candidate Institute for Computational and Engineering Sciences The University of Texas at Austin http://users.ices.utexas.edu/~vikram/ http://www.runforindia.org/runners/vikramg 
From: 蔡园武 <yuanwucai@gm...>  20120829 02:02:27

2012/8/29 Vikram Garg <simulationist@...>: > Hey Yuanwu, > Can you tell us what kind of linear solver options you > are using (the method, preconditioner options passed at runtime) ? I used the default option. (I don't know what's the default?) #EquationSystems n_systems()=1 System #0, "homogenization_elastic_3D_plate" Type "LinearImplicit" Variables="u" "v" "w" Finite Element Types="LAGRANGE" "LAGRANGE" "LAGRANGE" Approximation Orders="FIRST" "FIRST" "FIRST" n_dofs()=648 n_local_dofs()=648 n_constrained_dofs()=198 n_local_constrained_dofs()=198 n_vectors()=13 n_matrices()=1 DofMap Sparsity Average OnProcessor Bandwidth <= 84 Average OffProcessor Bandwidth <= 0 Maximum OnProcessor Bandwidth <= 144 Maximum OffProcessor Bandwidth <= 0 DofMap Constraints Number of DoF Constraints = 198 Average DoF Constraint Length= 1 > If you > are using the same stiffness matrix (K) and just changing the rhs (F), That's just what I'm doing. > then > you might actually be better of just using the LU preconditioner for one > load and reusing it for the others. How to set it in my code? Can you give me a simple example? Thanks for your answer! > > Thanks. > > On Tue, Aug 28, 2012 at 8:45 PM, 蔡园武 <yuanwucai@...> wrote: >> >> Hi, guys, >> I have a LinearImplicitSystem, solved with a sequential of different >> force vectors (rhs). >> Actually I defined different 'Fe' assemble function using a 'loadcase' >> indicator. In main function, I set the 'loadcase' value, call >> system.solve(), then es.reinit(), set a new 'loadcase' value, and >> solve() again. >> But I found that for some 'loadcase', the solver converged badly, like: >> >> loadcase1: Linear solver converged at step: 10736, final residual: >> 2.00331e21 >> loadcase2: Linear solver converged at step: 8549, final residual: >> 2.10685e21 >> loadcase3: Linear solver converged at step: 8, final residual: 1.38269e07 >> loadcase4: Linear solver converged at step: 0, final residual: 0.463112 >> loadcase5: Linear solver converged at step: 0, final residual: 0.463112 >> loadcase6: Linear solver converged at step: 0, final residual: 0.463112 >> >> In loadcase3, Linear solver converged at step 8, final resudual is not >> small enough. It's very strange that the loadcases after this didn't >> run? (converged at step 0?) The results are wrong and unbelievable. >> But if I solve loadcase3 after all the other loadcases, then loadcase >> 4,5,6 will be all right. Don't know why? >> >> Can I mannually control the final resudual tolerance? I did set a >> parameter in es: >> es.parameters.set<unsigned int> ("linear solver maximum >> iterations") = 20000; >> es.parameters.set<Real> ("linear solver tolerance") = TOLERANCE; >> >> Thanks for your help！ >>  >> Cai Yuanwu 蔡园武 >> Dept. of Engineering Mechanics, >> Dalian University of Technology, >> Dalian 116024, China >> >> >>  >> 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/ >> _______________________________________________ >> Libmeshusers mailing list >> Libmeshusers@... >> https://lists.sourceforge.net/lists/listinfo/libmeshusers > > > > >  > Vikram Garg > PhD Candidate > Institute for Computational and Engineering Sciences > The University of Texas at Austin > > http://users.ices.utexas.edu/~vikram/ > http://www.runforindia.org/runners/vikramg >  Cai Yuanwu 蔡园武 Dept. of Engineering Mechanics, Dalian University of Technology, Dalian 116024, China 
From: Vikram Garg <simulationist@gm...>  20120829 01:51:26

Hey Yuanwu, Can you tell us what kind of linear solver options you are using (the method, preconditioner options passed at runtime) ? If you are using the same stiffness matrix (K) and just changing the rhs (F), then you might actually be better of just using the LU preconditioner for one load and reusing it for the others. Thanks. On Tue, Aug 28, 2012 at 8:45 PM, 蔡园武 <yuanwucai@...> wrote: > Hi, guys, > I have a LinearImplicitSystem, solved with a sequential of different > force vectors (rhs). > Actually I defined different 'Fe' assemble function using a 'loadcase' > indicator. In main function, I set the 'loadcase' value, call > system.solve(), then es.reinit(), set a new 'loadcase' value, and > solve() again. > But I found that for some 'loadcase', the solver converged badly, like: > > loadcase1: Linear solver converged at step: 10736, final residual: > 2.00331e21 > loadcase2: Linear solver converged at step: 8549, final residual: > 2.10685e21 > loadcase3: Linear solver converged at step: 8, final residual: 1.38269e07 > loadcase4: Linear solver converged at step: 0, final residual: 0.463112 > loadcase5: Linear solver converged at step: 0, final residual: 0.463112 > loadcase6: Linear solver converged at step: 0, final residual: 0.463112 > > In loadcase3, Linear solver converged at step 8, final resudual is not > small enough. It's very strange that the loadcases after this didn't > run? (converged at step 0?) The results are wrong and unbelievable. > But if I solve loadcase3 after all the other loadcases, then loadcase > 4,5,6 will be all right. Don't know why? > > Can I mannually control the final resudual tolerance? I did set a > parameter in es: > es.parameters.set<unsigned int> ("linear solver maximum > iterations") = 20000; > es.parameters.set<Real> ("linear solver tolerance") = TOLERANCE; > > Thanks for your help！ >  > Cai Yuanwu 蔡园武 > Dept. of Engineering Mechanics, > Dalian University of Technology, > Dalian 116024, China > > >  > 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/ > _______________________________________________ > Libmeshusers mailing list > Libmeshusers@... > https://lists.sourceforge.net/lists/listinfo/libmeshusers >  Vikram Garg PhD Candidate Institute for Computational and Engineering Sciences The University of Texas at Austin http://users.ices.utexas.edu/~vikram/ http://www.runforindia.org/runners/vikramg 
From: 蔡园武 <yuanwucai@gm...>  20120829 01:45:21

Hi, guys, I have a LinearImplicitSystem, solved with a sequential of different force vectors (rhs). Actually I defined different 'Fe' assemble function using a 'loadcase' indicator. In main function, I set the 'loadcase' value, call system.solve(), then es.reinit(), set a new 'loadcase' value, and solve() again. But I found that for some 'loadcase', the solver converged badly, like: loadcase1: Linear solver converged at step: 10736, final residual: 2.00331e21 loadcase2: Linear solver converged at step: 8549, final residual: 2.10685e21 loadcase3: Linear solver converged at step: 8, final residual: 1.38269e07 loadcase4: Linear solver converged at step: 0, final residual: 0.463112 loadcase5: Linear solver converged at step: 0, final residual: 0.463112 loadcase6: Linear solver converged at step: 0, final residual: 0.463112 In loadcase3, Linear solver converged at step 8, final resudual is not small enough. It's very strange that the loadcases after this didn't run? (converged at step 0?) The results are wrong and unbelievable. But if I solve loadcase3 after all the other loadcases, then loadcase 4,5,6 will be all right. Don't know why? Can I mannually control the final resudual tolerance? I did set a parameter in es: es.parameters.set<unsigned int> ("linear solver maximum iterations") = 20000; es.parameters.set<Real> ("linear solver tolerance") = TOLERANCE; Thanks for your help！  Cai Yuanwu 蔡园武 Dept. of Engineering Mechanics, Dalian University of Technology, Dalian 116024, China 
From: Roy Stogner <roystgnr@ic...>  20120828 19:02:43

On Tue, 28 Aug 2012, David Knezevic wrote: > So we just leave points on a boundary that don't match up untouched? Works > for me. Right. find_neighbors() is designed to handle such as situation on meshes where the coarse mesh is conforming; if the coarse mesh is nonconforming then the user ought to be doing some fancy mortar method and not trying to stitch the topologies together anyway. Make sure you stick a libmesh_assert(is_serial) in here somewhere; getting a ParallelMesh union to work will probably require more complicated communications code than you want to bother with.  Roy 
From: David Knezevic <dknezevic@se...>  20120828 18:55:25

On 08/28/2012 02:52 PM, Roy Stogner wrote: > > On Tue, 28 Aug 2012, David Knezevic wrote: > >> OK, well I'd be happy to add the complete API, as long as we're >> happy to assume that the meshes are conforming on the interfaces, >> and that the user provides all the necessary connectivity info >> (essentially a "metamesh" data structure) > > Should it require a complicated data structure? User specifies two > boundary condition ids, library stitches together points from the one > that equal points from the other; user repeats with different ids as > necessary, then makes sure to prepare_for_use() before finishing. Yep, that's what I had in mind. The boundary IDs act as neighbor connectivity. So we just leave points on a boundary that don't match up untouched? Works for me. 
From: Roy Stogner <roystgnr@ic...>  20120828 18:52:17

On Tue, 28 Aug 2012, David Knezevic wrote: > OK, well I'd be happy to add the complete API, as long as we're > happy to assume that the meshes are conforming on the interfaces, > and that the user provides all the necessary connectivity info > (essentially a "metamesh" data structure) Should it require a complicated data structure? User specifies two boundary condition ids, library stitches together points from the one that equal points from the other; user repeats with different ids as necessary, then makes sure to prepare_for_use() before finishing.  Roy 
From: David Knezevic <dknezevic@se...>  20120828 18:48:19

On 08/28/2012 02:44 PM, Roy Stogner wrote: > > > On Tue, 28 Aug 2012, John Peterson wrote: > >> On Tue, Aug 28, 2012 at 12:16 PM, David Knezevic >> <dknezevic@...> wrote: >> >>> I agree, it'd be good to handle that in this function, but I was >>> originally >>> planning to do that in user code. The main reason for that is that I >>> wanted >>> to make use of extra info that I have, e.g. which mesh connects with >>> which >>> other mesh, and the boundary IDs of the mesh interfaces. I don't >>> really want >>> to put all the logic inside this union function... >> >> OK... I just worry that people will use this function and the >> resulting mesh without reading the documentation carefully, and get a >> nasty surprise. > > Likewise. I'd assumed that the MeshTools::_build_mesh_union() would > merely be a factoredapart building block for the more useful APIs. > If the more complete methods don't end up existing in library code (or > at least examples), I'm not sure how useful the bare building block > will be. OK, well I'd be happy to add the complete API, as long as we're happy to assume that the meshes are conforming on the interfaces, and that the user provides all the necessary connectivity info (essentially a "metamesh" data structure) 
From: David Knezevic <dknezevic@se...>  20120828 18:44:33

On 08/28/2012 02:39 PM, John Peterson wrote: > On Tue, Aug 28, 2012 at 12:16 PM, David Knezevic > <dknezevic@...> wrote: >> On 08/28/2012 01:23 PM, John Peterson wrote: >>> On Tue, Aug 28, 2012 at 10:34 AM, David Knezevic >>> <dknezevic@...> wrote: >>>> I'm doing some domain decomposition stuff, and I want to "stitch" >>>> distinct (conforming) meshes together into a single global mesh. The >>>> first step is to implement a "mesh union" operation, and I'll deal with >>>> merging "overlapping" nodes afterwards. >>>> >>>> To do this, I wrote a function build_mesh_union, which is pasted below. >>>> Worth adding to MeshTools, or somewhere like that? >>> I'd say definitely yes, but by "deal with merging overlapping nodes >>> afterwards" did you mean *this* function would eventually do that, or >>> that will be up to the user/another function? >>> >>> It seems to me that this function will be most useful if it can handle >>> the overlaps itself... >> >> I agree, it'd be good to handle that in this function, but I was originally >> planning to do that in user code. The main reason for that is that I wanted >> to make use of extra info that I have, e.g. which mesh connects with which >> other mesh, and the boundary IDs of the mesh interfaces. I don't really want >> to put all the logic inside this union function... > OK... I just worry that people will use this function and the > resulting mesh without reading the documentation carefully, and get a > nasty surprise. > I agree. I'll describe it in the header as: "Build a new mesh that is the union of all the meshes in \p meshes. This function currently does not attempt to merge the meshes, so in general the mesh that is returned contains meshes.size() disconnected submeshes." But I think it's also tricky to provide a general merging/stitching implementation. What if someone wants to take the union of nonconforming meshes, or even overlapping meshes? I think the current implementation at least has the advantage of simplicity, and users can then postprocess as necessary. 
From: Roy Stogner <roystgnr@ic...>  20120828 18:44:22

On Tue, 28 Aug 2012, John Peterson wrote: > On Tue, Aug 28, 2012 at 12:16 PM, David Knezevic > <dknezevic@...> wrote: > >> I agree, it'd be good to handle that in this function, but I was originally >> planning to do that in user code. The main reason for that is that I wanted >> to make use of extra info that I have, e.g. which mesh connects with which >> other mesh, and the boundary IDs of the mesh interfaces. I don't really want >> to put all the logic inside this union function... > > OK... I just worry that people will use this function and the > resulting mesh without reading the documentation carefully, and get a > nasty surprise. Likewise. I'd assumed that the MeshTools::_build_mesh_union() would merely be a factoredapart building block for the more useful APIs. If the more complete methods don't end up existing in library code (or at least examples), I'm not sure how useful the bare building block will be.  Roy 
From: John Peterson <jwpeterson@gm...>  20120828 18:39:55

On Tue, Aug 28, 2012 at 12:16 PM, David Knezevic <dknezevic@...> wrote: > > On 08/28/2012 01:23 PM, John Peterson wrote: >> >> On Tue, Aug 28, 2012 at 10:34 AM, David Knezevic >> <dknezevic@...> wrote: >>> >>> I'm doing some domain decomposition stuff, and I want to "stitch" >>> distinct (conforming) meshes together into a single global mesh. The >>> first step is to implement a "mesh union" operation, and I'll deal with >>> merging "overlapping" nodes afterwards. >>> >>> To do this, I wrote a function build_mesh_union, which is pasted below. >>> Worth adding to MeshTools, or somewhere like that? >> >> I'd say definitely yes, but by "deal with merging overlapping nodes >> afterwards" did you mean *this* function would eventually do that, or >> that will be up to the user/another function? >> >> It seems to me that this function will be most useful if it can handle >> the overlaps itself... > > > I agree, it'd be good to handle that in this function, but I was originally > planning to do that in user code. The main reason for that is that I wanted > to make use of extra info that I have, e.g. which mesh connects with which > other mesh, and the boundary IDs of the mesh interfaces. I don't really want > to put all the logic inside this union function... OK... I just worry that people will use this function and the resulting mesh without reading the documentation carefully, and get a nasty surprise.  John 
From: David Knezevic <dknezevic@se...>  20120828 18:16:41

On 08/28/2012 01:23 PM, John Peterson wrote: > On Tue, Aug 28, 2012 at 10:34 AM, David Knezevic > <dknezevic@...> wrote: >> I'm doing some domain decomposition stuff, and I want to "stitch" >> distinct (conforming) meshes together into a single global mesh. The >> first step is to implement a "mesh union" operation, and I'll deal with >> merging "overlapping" nodes afterwards. >> >> To do this, I wrote a function build_mesh_union, which is pasted below. >> Worth adding to MeshTools, or somewhere like that? > I'd say definitely yes, but by "deal with merging overlapping nodes > afterwards" did you mean *this* function would eventually do that, or > that will be up to the user/another function? > > It seems to me that this function will be most useful if it can handle > the overlaps itself... I agree, it'd be good to handle that in this function, but I was originally planning to do that in user code. The main reason for that is that I wanted to make use of extra info that I have, e.g. which mesh connects with which other mesh, and the boundary IDs of the mesh interfaces. I don't really want to put all the logic inside this union function... David 
From: David Knezevic <dknezevic@se...>  20120828 18:11:03

On 08/28/2012 01:21 PM, Roy Stogner wrote: > > On Tue, 28 Aug 2012, David Knezevic wrote: > >> I'm doing some domain decomposition stuff, and I want to "stitch" >> distinct (conforming) meshes together into a single global mesh. The >> first step is to implement a "mesh union" operation, and I'll deal with >> merging "overlapping" nodes afterwards. >> >> To do this, I wrote a function build_mesh_union, which is pasted below. >> Worth adding to MeshTools, or somewhere like that? > > In my opinion, definitely worth adding to MeshTools. OK, great. > Any chance it could be defined to work with ranges and output > iterators rather than just MeshBases? I'm pretty proud of the way > MeshCommunication::allgather() works now, for example. You'd then > have the ability to do things like stitching together selected > subdomains "for free". Sounds good to me, but I'd have to look at some example code in the library to figure out how to implement that. How about I check in what I have now, and then we add the ranges and output iterators later? David 
From: John Peterson <jwpeterson@gm...>  20120828 17:24:21

On Tue, Aug 28, 2012 at 10:34 AM, David Knezevic <dknezevic@...> wrote: > I'm doing some domain decomposition stuff, and I want to "stitch" > distinct (conforming) meshes together into a single global mesh. The > first step is to implement a "mesh union" operation, and I'll deal with > merging "overlapping" nodes afterwards. > > To do this, I wrote a function build_mesh_union, which is pasted below. > Worth adding to MeshTools, or somewhere like that? I'd say definitely yes, but by "deal with merging overlapping nodes afterwards" did you mean *this* function would eventually do that, or that will be up to the user/another function? It seems to me that this function will be most useful if it can handle the overlaps itself...  John 