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}
(100) 
_{Nov}

_{Dec}

S  M  T  W  T  F  S 







1

2

3
(4) 
4
(1) 
5
(1) 
6

7
(4) 
8
(1) 
9

10

11

12
(1) 
13

14

15

16

17

18
(4) 
19
(3) 
20

21

22
(4) 
23
(6) 
24

25

26
(5) 
27
(2) 
28
(6) 
29
(3) 
30
(1) 
31
(1) 





From: <tim@ce...>  20060731 06:13:38

Thanks, Derek, for your reply. Could please someone who has write access to the CVS repository check your changes in, so that other peeople don't run into the same problems? On Sat, 29 Jul 2006, Derek Gaston wrote: >> 2.) The constructor that takes "unsigend int var" (rather than >> "std::vector<unsigend int> vars") initializes the member _system_vars >> with "std::vector<unsigned int>(var)" which is a vector of var >> components, all being zero, rather than "std::vector<unsigned >> int>(1,var)". My workaround is to use the other constructor. > > I don't know about this... I must have used the other constructor all > along. This is probably something to be looked into. It's a very easy task: In line 65 of mesh_function.C, the "(var)" has to be changed into "(1,var)". Could someone please also check in this change? > Also... to everyone else.... why is it that mesh_funciton.C is in > numerics when mesh_function.h is in mesh? I personally vote for both > being in mesh! I agree. Best Regards, Tim 
From: Steffen Petersen <steffen.petersen@tu...>  20060730 08:44:23

David Xu schrieb: > Hi All, > > I'm compiled libMesh with Slepc/ARPACK for generalized real eigenvalue > problems. I was wondering if it's possible to a guess eigenvalue and ask > ARPACK to search for solutions near or greater to the guess? The current You will have to apply a spectral transformation in order to shift the spectrum to the frequency range of interest. The easiest way to apply the shift within your libMesh program is to make use of the slepc command line option "st_shift". In chapter 3 of the slepc manual you will find more on this issue. Steffen > options I know of are: > > LARGEST_MAGNITUDE, SMALLEST_MAGNITUDE, LARGEST_REAL, SMALLEST_REAL, > LARGEST_IMAGINARY, SMALLEST_IMAGINARY > > Thanks, > > David > > >  > >  > Take Surveys. Earn Cash. Influence the Future of IT > Join SourceForge.net's Techsay panel and you'll get the chance to share your > opinions on IT & business topics through brief surveys  and earn cash > http://www.techsay.com/default.php?page=join.php&p=sourceforge&CID=DEVDEV > > >  > > _______________________________________________ > Libmeshusers mailing list > Libmeshusers@... > https://lists.sourceforge.net/lists/listinfo/libmeshusers 
From: Derek Gaston <friedmud@gm...>  20060729 19:37:45

> 1.) The MeshFunction class cannot be instantiated since the clear() > function is abstract. My workaround: I use a derived class > MyMeshFunction in which I implement a clear() function that does > nothing. This is essentially what I did... except I just added it to the actual class itself... in the header file. I've attached the diff. > > 2.) The constructor that takes "unsigend int var" (rather than > "std::vector<unsigend int> vars") initializes the member _system_vars > with "std::vector<unsigned int>(var)" which is a vector of var > components, all being zero, rather than "std::vector<unsigned > int>(1,var)". My workaround is to use the other constructor. I don't know about this... I must have used the other constructor all along. This is probably something to be looked into. > > 3.) When calling the init() function, I get the message "ERROR: > Already initialized! Will ignore this call...", although I definitely > call this function only once. All you have to do is comment out this>_point_locator>init(); in mesh_function.C. It looks like the point_locator class does this automatically now... so I think this is crusty code that wasn't updated when point_locator was... the problem is it is trying to initialize it twice... I've attached my diff for that as well. > Best Regards, > > Tim Hope that helps! Also... to everyone else.... why is it that mesh_funciton.C is in numerics when mesh_function.h is in mesh? I personally vote for both being in mesh! Derek 
From: Roy Stogner <roystgnr@ic...>  20060729 08:54:52

On Fri, 28 Jul 2006, David Xu wrote: > I was wondering how to write down an equation containing mixed derivative in > libMesh's element matrix (Ke) building step? > > For exmaple, Uxx+Uyy+Uxy (U is 2nd order derivative) If you need mixed second derivatives of a shape function, take a look at example 15, which solves a biharmonic equation by a direct Galerkin method which computes laplacians at each quadrature point. If you need derivatives of your solution, take a look at how solutions and their derivatives are accumulated in example 13  you'd do the same thing but with Tensor/RealTensor classes instead of Gradient/RealGradient.  Roy 
From: Derek Gaston <friedmud@gm...>  20060729 06:50:30

Tim, I am currently using MeshFunction... and I did run into the problems you spoke about... and I've worked around each in similar ways. It's getting pretty late now... so I will respond tomorrow with more details, but I just wanted to get a response out there. To everyone else: I was going to report these problems at some point, but at the time I was knee deep in my own problems and was just trying to get stuff working... hopefully after my comments tomorrow we can get these small problems resolved. Derek On 7/28/06, Tim Kr=F6ger <tim@...> wrote: > Dear all, > > During my work with libMesh, I encountered the following three > problems with the MeshFunction class: > > 1.) The MeshFunction class cannot be instantiated since the clear() > function is abstract. My workaround: I use a derived class > MyMeshFunction in which I implement a clear() function that does > nothing. > > 2.) The constructor that takes "unsigend int var" (rather than > "std::vector<unsigend int> vars") initializes the member _system_vars > with "std::vector<unsigned int>(var)" which is a vector of var > components, all being zero, rather than "std::vector<unsigned > int>(1,var)". My workaround is to use the other constructor. > > 3.) When calling the init() function, I get the message "ERROR: > Already initialized! Will ignore this call...", although I definitely > call this function only once. > > Can anyone comment on this? > > Best Regards, > > Tim > >  > Take Surveys. Earn Cash. Influence the Future of IT > Join SourceForge.net's Techsay panel and you'll get the chance to share y= our > opinions on IT & business topics through brief surveys  and earn cash > http://www.techsay.com/default.php?page=3Djoin.php&p=3Dsourceforge&CID=3D= DEVDEV > _______________________________________________ > Libmeshusers mailing list > Libmeshusers@... > https://lists.sourceforge.net/lists/listinfo/libmeshusers > 
From: Roy Stogner <roystgnr@ic...>  20060728 22:19:58

This message should only apply to people using: BERNSTEIN elements with p > 2 CLOUGH elements HIERARCHIC elements with p > 2 SZABAB elements with p > 2 Currently, we choose an arbitrary orientation for certain basis functions by examining the relative global node numbering of the nodes on the appropriate element edge or face. I would like to change that procedure, so that instead we choose basis orientations based on the relative locations of those nodes. The advantage is that this will enable more optimizations in the FE class (cutting as much as 30% off the runtime of one of my application codes). The disadvantage is that if anyone has any saved solution .xda files using these elements, those solutions will become corrupted in future libMesh versions. If anyone has such saved solution files that they can't bear to part with or recalculate, email me and I'll try and work out a way to get them converted. If not (and oh, I'm hoping there's not), then just smile and look forward to your code running a little faster with libMesh 0.6.0.  Roy Stogner 
From: David Xu <dxu@my...>  20060728 22:17:09

Hi All, I was wondering how to write down an equation containing mixed derivative in libMesh's element matrix (Ke) building step? For exmaple, Uxx+Uyy+Uxy (U is 2nd order derivative) Thanks! David 
From: Roy Stogner <roystgnr@ic...>  20060728 21:48:24

On Wed, 26 Jul 2006, David Xu wrote: > On 7/26/06, Roy Stogner <roystgnr@...> wrote: > >> Loop over all the elements. On each geometric element, reinit your FE >> object with a grid quadrature rule whose order at least equals your >> polynomial degree. Get the XYZ coordinates and the solution value at >> each quadrature point, and output them. > > I do this after I load the eigenvector solutions back to libMesh, correct? Yes. > If I reinit FE object with a grid quadrature rule, how can Imake sure the > XYZ coordinates at each quadrature point are the ones correspondent to the > loaded eigenvectors? You don't. You're not plotting eigenvector coefficients, you're plotting eigenfunction points. Depending on the finite element space you use, your coefficients may correspond to function points, but they may also correspond to mixed second derivatives, projection magnitudes onto orthogonal function spaces, or the phase of the moon! So forget about trying to plot the coefficients themselves, and just worry about how to plot the function they represent. > Can I just ouput the xyz at each quadrature point at the element matrices > (Ke, Me) assembly step? Output them to a file and concatenate with the > eigenvector solutions? Will that work? No. >> That will output vertex and edge points multiple times, but the plot >> should look the same. If your elements are so large that the plot >> still looks faceted, use a quadrature rule with more points. > > Same question here, if I use a quadrature rule different from the one used > to assemble Ke, Me, will the xyz coordinates at the each quadrature point be > correctlycorrespondent to the eigenvector solutions with the size of global > DOFs? No, they won't. That's why you don't just get the xyz coordinates at each quadrature point, you get the solution value at each quadrature point too.  Roy 
From: David Xu <dxu@my...>  20060728 20:54:01

Hi All, I'm compiled libMesh with Slepc/ARPACK for generalized real eigenvalue problems. I was wondering if it's possible to a guess eigenvalue and ask ARPACK to search for solutions near or greater to the guess? The current options I know of are: LARGEST_MAGNITUDE, SMALLEST_MAGNITUDE, LARGEST_REAL, SMALLEST_REAL, LARGEST_IMAGINARY, SMALLEST_IMAGINARY Thanks, David 
From: <tim@ce...>  20060728 06:33:07

Dear all, During my work with libMesh, I encountered the following three problems with the MeshFunction class: 1.) The MeshFunction class cannot be instantiated since the clear() function is abstract. My workaround: I use a derived class MyMeshFunction in which I implement a clear() function that does nothing. 2.) The constructor that takes "unsigend int var" (rather than "std::vector<unsigend int> vars") initializes the member _system_vars with "std::vector<unsigned int>(var)" which is a vector of var components, all being zero, rather than "std::vector<unsigned int>(1,var)". My workaround is to use the other constructor. 3.) When calling the init() function, I get the message "ERROR: Already initialized! Will ignore this call...", although I definitely call this function only once. Can anyone comment on this? Best Regards, Tim 
From: Roy Stogner <roystgnr@ic...>  20060728 03:47:30

On Wed, 26 Jul 2006, David Xu wrote: > Thanks. The reason why I need xyz of the DOFs is that my system matrices are > solved by an eigenvalue solver outside libmesh and the size of the > eigenvector solution produced by the solver equals to the dimension of the > system matrices (or DOFs). I'd like to plot the eigenvector solutions, thus > I need to the correspondent xyz coordinates of each of the eigenvector value > for visualization purpose. No, you don't. Those XYZ coordinates *overlap*  each Hermite point will have 2 degrees of freedom in 1D, 4 in 2D, and 8 in 3D. You can throw away all but the first degree of freedom on each node and plot each piecewise cubic element as a piecewise linear, but if you want an accurate representation of your results you'll either need plotting software that can handle higher polynomials or you'll need to subdivide your mesh elements. The most general way you can do things is to load your eigenvector solutions back into libMesh (preferably into the same running process that gave you the original matrix, so there's no question of node renumbering). Then you can use the libMesh plotting functions (assuming you're happy with the limitations of those output formats) or use the libMesh FE objects to get an arbitrarily dense cloud of points for your own plotting software.  Roy Stogner 
From: Roy Stogner <roystgnr@ic...>  20060727 14:53:43

On Wed, 26 Jul 2006, David Xu wrote: > I kinda get it, but still vague about how to implement it. Once I load the > eigenfunction value points from file to libMesh using > NumericVector<>::set(), how does libMesh provide me with the eigenfunction > point at requested quadrature point? I think I'm missing something here. I > know how to get xyz at requested quadrature point from the examples. See example 13 for how to get solution values at quadrature points. You just sum the basis function falues weighted by the local coefficients from the solution.  Roy 
From: David Xu <dxu@my...>  20060727 03:31:28

On 7/26/06, Roy Stogner <roystgnr@...> wrote: > > On Wed, 26 Jul 2006, David Xu wrote: > > > On 7/26/06, Roy Stogner <roystgnr@...> wrote: > > > >> Loop over all the elements. On each geometric element, reinit your FE > >> object with a grid quadrature rule whose order at least equals your > >> polynomial degree. Get the XYZ coordinates and the solution value at > >> each quadrature point, and output them. I kinda get it, but still vague about how to implement it. Once I load the eigenfunction value points from file to libMesh using NumericVector<>::set(), how does libMesh provide me with the eigenfunction point at requested quadrature point? I think I'm missing something here. I know how to get xyz at requested quadrature point from the examples. > > > I do this after I load the eigenvector solutions back to libMesh, > correct? > > Yes. > > > If I reinit FE object with a grid quadrature rule, how can Imake sure > the > > XYZ coordinates at each quadrature point are the ones correspondent to > the > > loaded eigenvectors? > > You don't. You're not plotting eigenvector coefficients, you're > plotting eigenfunction points. Depending on the finite element space > you use, your coefficients may correspond to function points, but they > may also correspond to mixed second derivatives, projection > magnitudes onto orthogonal function spaces, or the phase of the moon! > So forget about trying to plot the coefficients themselves, and just > worry about how to plot the function they represent. > > > Can I just ouput the xyz at each quadrature point at the element > matrices > > (Ke, Me) assembly step? Output them to a file and concatenate with the > > eigenvector solutions? Will that work? > > No. > > >> That will output vertex and edge points multiple times, but the plot > >> should look the same. If your elements are so large that the plot > >> still looks faceted, use a quadrature rule with more points. > > > > Same question here, if I use a quadrature rule different from the one > used > > to assemble Ke, Me, will the xyz coordinates at the each quadrature > point be > > correctlycorrespondent to the eigenvector solutions with the size of > global > > DOFs? > > No, they won't. That's why you don't just get the xyz coordinates at > each quadrature point, you get the solution value at each quadrature > point too. >  > Roy > 
From: David Xu <dxu@my...>  20060726 20:34:59

On 7/26/06, Roy Stogner <roystgnr@...> wrote: > > On Wed, 26 Jul 2006, David Xu wrote: > > > On 7/26/06, Roy Stogner <roystgnr@...> wrote: > > > > I see. Do you any of the plotting software that can handle higher > > polynomials? > > No. If I knew of any I'd be using it myself. ;) As it is I usually > get good enough plots from the bisected piecewise linear plots that > libMesh's output functions produce by default. > > > I remember you mentioned lagrange only support up to 2nd order. If > that's > > correct, > > That's correct. > > > I may also want to try higher order hierarchic (up to p=5?). > > The polynomial degrees supported by each finite element class depend > on the geometric element. Check out the definition of > FEInterface::max_order() if you want to see what each supports. > > On HEX27 elements, the hierarchics support "unlimited" polynomial > degree  but ill conditioning and floating point error will probably > screw up your solution around p=11. > > > Will hierarchic (p>2) have the same node renumbering problem for > > getting xyz? > > Wait wait wait  you're talking about two problems here. When I > mentioned node renumbering before, I was talking about the difficulty > of loading old solutions into new meshes. The problem in plotting > solutions on nonlagrange elements is different, and yes, the > hierarchics will be equally tricky. > > >> The most general way you can do things is to load your eigenvector > >> solutions back into libMesh (preferably into the same running process > >> that gave you the original matrix, so there's no question of node > >> renumbering). Then you can use the libMesh plotting functions > >> (assuming you're happy with the limitations of those output formats) > >> or use the libMesh FE objects to get an arbitrarily dense cloud of > >> points for your own plotting software. > > > > That sounds like a project to me. Would you please point me to the > possible > > classes/functions I need to load the eigenvector solutions back into > libMesh > > Yes: NumericVector<>::set() We don't have any "load numeric vector > from file" function, but it shouldn't be too hard to read in a file > and set the coefficients yourself. > > > and to plot them? Currently, for the 2nd order lagrange, I just > concatenate > > the node xyz locations with the the eigenvector solutions, and then plot > > > them in Matlab. > > Okay. If Matlab doesn't have any particular requirements for the > ordering of the points you plot, then it's easy: No, the ordering of the points is not important in matlab plotting, I think. Basically I loaded a matrix that contains x, y, z, eigenvectors as inividual columns and used surface/mesh/line plot. The only requirement is the each xyz physical location is correspondent to the eigenvector solution. Loop over all the elements. On each geometric element, reinit your FE > object with a grid quadrature rule whose order at least equals your > polynomial degree. Get the XYZ coordinates and the solution value at > each quadrature point, and output them. I do this after I load the eigenvector solutions back to libMesh, correct? If I reinit FE object with a grid quadrature rule, how can Imake sure the XYZ coordinates at each quadrature point are the ones correspondent to the loaded eigenvectors? Can I just ouput the xyz at each quadrature point at the element matrices (Ke, Me) assembly step? Output them to a file and concatenate with the eigenvector solutions? Will that work? That will output vertex and edge points multiple times, but the plot > should look the same. If your elements are so large that the plot > still looks faceted, use a quadrature rule with more points. Same question here, if I use a quadrature rule different from the one used to assemble Ke, Me, will the xyz coordinates at the each quadrature point be correctlycorrespondent to the eigenvector solutions with the size of global DOFs? Thanks, David 
From: Roy Stogner <roystgnr@ic...>  20060726 17:59:10

On Wed, 26 Jul 2006, David Xu wrote: > On 7/26/06, Roy Stogner <roystgnr@...> wrote: > > I see. Do you any of the plotting software that can handle higher > polynomials? No. If I knew of any I'd be using it myself. ;) As it is I usually get good enough plots from the bisected piecewise linear plots that libMesh's output functions produce by default. > I remember you mentioned lagrange only support up to 2nd order. If that's > correct, That's correct. > I may also want to try higher order hierarchic (up to p=5?). The polynomial degrees supported by each finite element class depend on the geometric element. Check out the definition of FEInterface::max_order() if you want to see what each supports. On HEX27 elements, the hierarchics support "unlimited" polynomial degree  but ill conditioning and floating point error will probably screw up your solution around p=11. > Will hierarchic (p>2) have the same node renumbering problem for > getting xyz? Wait wait wait  you're talking about two problems here. When I mentioned node renumbering before, I was talking about the difficulty of loading old solutions into new meshes. The problem in plotting solutions on nonlagrange elements is different, and yes, the hierarchics will be equally tricky. >> The most general way you can do things is to load your eigenvector >> solutions back into libMesh (preferably into the same running process >> that gave you the original matrix, so there's no question of node >> renumbering). Then you can use the libMesh plotting functions >> (assuming you're happy with the limitations of those output formats) >> or use the libMesh FE objects to get an arbitrarily dense cloud of >> points for your own plotting software. > > That sounds like a project to me. Would you please point me to the possible > classes/functions I need to load the eigenvector solutions back into libMesh Yes: NumericVector<>::set() We don't have any "load numeric vector from file" function, but it shouldn't be too hard to read in a file and set the coefficients yourself. > and to plot them? Currently, for the 2nd order lagrange, I just concatenate > the node xyz locations with the the eigenvector solutions, and then plot > them in Matlab. Okay. If Matlab doesn't have any particular requirements for the ordering of the points you plot, then it's easy: Loop over all the elements. On each geometric element, reinit your FE object with a grid quadrature rule whose order at least equals your polynomial degree. Get the XYZ coordinates and the solution value at each quadrature point, and output them. That will output vertex and edge points multiple times, but the plot should look the same. If your elements are so large that the plot still looks faceted, use a quadrature rule with more points.  Roy 
From: Roy Stogner <roystgnr@ic...>  20060726 17:13:39

On Wed, 26 Jul 2006, David Xu wrote: > After I used 3rd order hermite shape function, I realized the resulting > dimension of the assembled system matrices is not equal to the number of > mesh nodes as in the case of using 2nd order Lagrange. The dimension equals > the number of DOFs. I was wondering how to get xyz physical location of each > of the global DOF points. Is there an iterator available like the one for > mesh nodes? No, there isn't. If you want to be dangerous, you can use the DOFObject interface of a node in the mesh of an initialized EquationSystems to view every global degree of freedom index on that node. Keep in mind that for many elements (including the HERMITE elements for p=4 and above, I think) there are also degrees of freedom associated with the Elem itself. Why do you need the physical location of nonLagrange degrees of freedom? That's kind of unusual.  Roy Stogner 
From: David Xu <dxu@my...>  20060726 15:52:05

Roy, Thanks. The reason why I need xyz of the DOFs is that my system matrices are solved by an eigenvalue solver outside libmesh and the size of the eigenvector solution produced by the solver equals to the dimension of the system matrices (or DOFs). I'd like to plot the eigenvector solutions, thus I need to the correspondent xyz coordinates of each of the eigenvector value for visualization purpose. It was easy to do for the 2nd order Lagrange because the size of the eigenvector = n_nodes, now I'm stuck at plotting the results from 3rd order hermite, which have much better accuracy. David On 7/26/06, Roy Stogner <roystgnr@...> wrote: > > On Wed, 26 Jul 2006, David Xu wrote: > > > After I used 3rd order hermite shape function, I realized the resulting > > dimension of the assembled system matrices is not equal to the number of > > mesh nodes as in the case of using 2nd order Lagrange. The dimension > equals > > the number of DOFs. I was wondering how to get xyz physical location of > each > > of the global DOF points. Is there an iterator available like the one > for > > mesh nodes? > > No, there isn't. If you want to be dangerous, you can use the > DOFObject interface of a node in the mesh of an initialized > EquationSystems to view every global degree of freedom index on that > node. Keep in mind that for many elements (including the HERMITE > elements for p=4 and above, I think) there are also degrees of freedom > associated with the Elem itself. > > Why do you need the physical location of nonLagrange degrees of > freedom? That's kind of unusual. >  > Roy Stogner > 
From: David Xu <dxu@my...>  20060726 10:53:30

Hi All, After I used 3rd order hermite shape function, I realized the resulting dimension of the assembled system matrices is not equal to the number of mesh nodes as in the case of using 2nd order Lagrange. The dimension equals the number of DOFs. I was wondering how to get xyz physical location of each of the global DOF points. Is there an iterator available like the one for mesh nodes? Thanks, David 
From: Roy Stogner <roystgnr@ic...>  20060723 14:06:29

On Sun, 23 Jul 2006, David Knezevic wrote: > My main stumbling block is how to implement this loading > and storing of solution vectors in parallel. Yeah, mine too. ;) > It looks to me like the function "update_global_solution" in System > would enable me to store the solution computed on multiple > processors in a localized vector, which solves the "storing" part of > my problem. However, I can't see an easy way of broadcasting a > localized solution vector to the various processors. Is there > functionality for this kind of thing in libMesh, or would it require > specialized MPI code? Well, System::reinit() updates the distributed solution vector from the localized current_local_solution vector. I don't see any standard functionality for doing that with a different localized vector, but you could either use current_local_solution (that's what it's there for) or rewrite the relevant reinit() code (it's short and simple).  Roy 
From: David Knezevic <david.knezevic@ba...>  20060723 13:35:23

Hi all, I'm solving a PDE using an ADIlike method (which is possible because the PDE is amenable to operatorsplitting) in which I store the solution of the PDE in a matrix and sequentially pick out either a row or column of the matrix, load this vector into a libMesh equation system, update it and then store it again in the matrix. This reduces the problem I'm solving from a single 6dimensional problem to a series of 3D problems, which is a significant gain. Anyway, I have this working well in a singleprocessor libMesh implementation, but as you can imagine it's a very computationally intensive exercise so I'd like to move to implementing this in parallel. My main stumbling block is how to implement this loading and storing of solution vectors in parallel. It looks to me like the function "update_global_solution" in System would enable me to store the solution computed on multiple processors in a localized vector, which solves the "storing" part of my problem. However, I can't see an easy way of broadcasting a localized solution vector to the various processors. Is there functionality for this kind of thing in libMesh, or would it require specialized MPI code? Thanks for the help. Regards, David 
From: Roy Stogner <roystgnr@ic...>  20060723 03:02:49

On Sat, 22 Jul 2006, David Xu wrote: > On 7/22/06, Roy Stogner <roystgnr@...> wrote: >> >> On Sat, 22 Jul 2006, David Xu wrote: >> >> > Are you saying HEX27 and HEX8 will give same level of accuracy in the >> > solutions if I don't need quadratic mapping functions? >> >> Yes. > > Just ccurious. Does the same rule apply to other types of element? What > about Tet, Tri, Quad and Prism. So, is the level of solution accuracy > independent with the number of node within the same type of element? The level of solution accuracy is independent of the number of geometric nodes... but libMesh reuses geometric nodes to store degrees of the freedom that have the same topological connectivity, so you usually still need secondorder nodes even if you aren't fitting a secondorder geometry. As far as I know, the HERMITE elements and the two discontinuous elements are the only way to get better than linear approximations on linear geometric elements. If you try to use finite elements on geometric elements that don't support them, however, you won't just get reduced accuracy, your code will exit with an error. > What about the difference between different element types in terms > of the effect on the quality of the solutions? You can get better solutions (better conditioned matrices, at least) from quadratic elements if you use a mesh smoother that takes advantage of them. Mostly, though, you only need higher order geometric elements to better fit curved domain boundaries. >> Almost certainly it is. None of our finite element classes are as >> optimized as they should be, but I think the Hermite elements may be >> worse than average. >> >> Keep in mind, too, that even if they were equally optimized, the >> Hermite assembly would be more expensive. If you're using the default >> quadrature order (which is designed for nonlinear problems and may be >> gross overkill for you) then I think quadratic hexes will be >> calculating at 27 points and cubic hexes will be calculating at 64. > > That explains why the ouput filesize from hermite is much larger than > lagrange. No, it doesn't. I'm talking about quadrature points here, and the size of your final matrix is (with few exceptions) independent of the quadrature rule you use to calculate it. I can see why that's confusing, though: by coincidence the number of quadrature points is the same as the number of local DoFs for both elements here. Of course, it's not just the quadrature rule that's important. Having 64 local DoFs instead of 27 also increases calculation time. Finally, on uniform meshes Hermite cube DoFs usually couple to 216 DoFs rather than 27, 45, or 125, which is probably what's increasing your output file size. > Does that mean, even the matrix dimension is the same, but hermit > produces more entries in the matrix, thus it's less sparse? Yes. Increasing polynomial order requires more bandwidth, so does increaing continuity, and going from quadratic Lagrange to cubic Hermite does both at once.  Roy 
From: David Xu <dxu@my...>  20060723 02:31:40

On 7/22/06, Roy Stogner <roystgnr@...> wrote: > > On Sat, 22 Jul 2006, David Xu wrote: > > > On 7/22/06, Roy Stogner <roystgnr@...> wrote: > > > > Are you saying HEX27 and HEX8 will give same level of accuracy in the > > solutions if I don't need quadratic mapping functions? > > Yes. Just ccurious. Does the same rule apply to other types of element? What about Tet, Tri, Quad and Prism. So, is the level of solution accuracy independent with the number of node within the same type of element? What about the difference between different element types in terms of the effect on the quality of the solutions? > It's for solving 2 system matrices in a realvalued generalized eigenvalue > > problem. > > Okay, then the simplest workaround is clear: only build one matrix at > a time. As long as you're just writing them both out to files > anyway, there's no reason you need them both in RAM simultaneously. Yes, that's a great idea. > I went back and tried (30x30x30, HEX27, 3rd order HERMITE). This > > time it didn't blew out the memory, but did take significantly > > longer time to assemble the matrices than (30x30x30, HEX27, 2nd > > order LANGRANGE). Maybe HERMITE is the problem? > > Almost certainly it is. None of our finite element classes are as > optimized as they should be, but I think the Hermite elements may be > worse than average. > > Keep in mind, too, that even if they were equally optimized, the > Hermite assembly would be more expensive. If you're using the default > quadrature order (which is designed for nonlinear problems and may be > gross overkill for you) then I think quadratic hexes will be > calculating at 27 points and cubic hexes will be calculating at 64. That explains why the ouput filesize from hermite is much larger than lagrange. Does that mean, even the matrix dimension is the same, but hermit produces more entries in the matrix, thus it's less sparse? > How to assemble the element matrices outside of libmesh? > > Basically you'd almost do what libMesh does: create a big empty matrix > of the appropriate sparsity pattern, then loop though all the element > matrices and add their entries after looking up the global index for > each local degree of freedom. > > The only difference is that because you know you've got a uniform > grid, you could do that local>global lookup with a few equations > instead of the big data structures that general unstructured grids > require. > > > Do you know any existing code/program that can do that? So, that > > would be, output each element matrix to a file and a problem should > > be able to read in all the element matrices from the file and > > assemble them to system matrices. I'm definitely interested if this > > is doable. > > It's definitely doable, but I don't know of any existing code to do > it. I could script it up in Matlab pretty easily, but the Matlab > sparse matrix format sucks and so I wouldn't want to work with the > result. I might try it using python/numpy/scipy. Thanks for the great tips! David 
From: Roy Stogner <roystgnr@ic...>  20060723 01:22:43

On Sat, 22 Jul 2006, David Xu wrote: > On 7/22/06, Roy Stogner <roystgnr@...> wrote: > > Are you saying HEX27 and HEX8 will give same level of accuracy in the > solutions if I don't need quadratic mapping functions? Yes. > It's for solving 2 system matrices in a realvalued generalized eigenvalue > problem. Okay, then the simplest workaround is clear: only build one matrix at a time. As long as you're just writing them both out to files anyway, there's no reason you need them both in RAM simultaneously. > I went back and tried (30x30x30, HEX27, 3rd order HERMITE). This > time it didn't blew out the memory, but did take significantly > longer time to assemble the matrices than (30x30x30, HEX27, 2nd > order LANGRANGE). Maybe HERMITE is the problem? Almost certainly it is. None of our finite element classes are as optimized as they should be, but I think the Hermite elements may be worse than average. Keep in mind, too, that even if they were equally optimized, the Hermite assembly would be more expensive. If you're using the default quadrature order (which is designed for nonlinear problems and may be gross overkill for you) then I think quadratic hexes will be calculating at 27 points and cubic hexes will be calculating at 64. > How to assemble the element matrices outside of libmesh? Basically you'd almost do what libMesh does: create a big empty matrix of the appropriate sparsity pattern, then loop though all the element matrices and add their entries after looking up the global index for each local degree of freedom. The only difference is that because you know you've got a uniform grid, you could do that local>global lookup with a few equations instead of the big data structures that general unstructured grids require. > Do you know any existing code/program that can do that? So, that > would be, output each element matrix to a file and a problem should > be able to read in all the element matrices from the file and > assemble them to system matrices. I'm definitely interested if this > is doable. It's definitely doable, but I don't know of any existing code to do it. I could script it up in Matlab pretty easily, but the Matlab sparse matrix format sucks and so I wouldn't want to work with the result.  Roy 
From: David Xu <dxu@my...>  20060723 00:15:52

On 7/22/06, Roy Stogner <roystgnr@...> wrote: > > On Sat, 22 Jul 2006, David Xu wrote: > > > I was trying to assemble the stiffness and mass matrices on a dense mesh > > with 531441 nodes (40x40x40, HEX27, 3rd, HERMITE) > > Can I suggest trying HEX8? The HERMITE elements are unique among our > higher order elements in that all their degrees of freedom are > topologically associated with mesh vertices, so unless you need > quadratic mapping functions you don't need HEX27 elements. Those > nodes aren't reponsible for most of your memory use (the system matrix > is), but every megabyte helps. Are you saying HEX27 and HEX8 will give same level of accuracy in the solutions if I don't need quadratic mapping functions? > and I ran into "Out of memory" problem at the line of: > > > > equation_systems.init(); > > > > > > Here's the error message: > > > > [0]PETSC ERROR: PetscMallocAlign() line 62 in src/sys/memory/mal.c > > [0]PETSC ERROR: Out of memory. This could be due to allocating > > [0]PETSC ERROR: too large an object or bleeding by not properly > > [0]PETSC ERROR: destroying unneeded objects. > > [0]PETSC ERROR: Memory allocated 1380430780 Memory used by process > 886095872 > > [0]PETSC ERROR: Try running with malloc_dump or malloc_log for info. > > [0]PETSC ERROR: Memory requested 907039528! > > [0]PETSC ERROR: PetscTrMallocDefault() line 191 in src/sys/memory/mtr.c > > [0]PETSC ERROR: MatSeqAIJSetPreallocation_SeqAIJ() line 2735 in > > src/mat/impls/aij/seq/aij.c > > [0]PETSC ERROR: MatCreateSeqAIJ() line 2621 in > src/mat/impls/aij/seq/aij.c > > [0]PETSC ERROR: User provided function() line 137 in > > unknowndirectory/src/numerics/petsc_matrix.C > > [unset]: aborting job: > > application called MPI_Abort(comm=0x84000000, 1)  process 0 > > > > My question is: is it possible to assemble the matrices without having > to > > initilize the equation system? > > I'm afraid not. You could initialize a finite element object and > evaluate the element matrices, but putting them into the system matrix > requires that matrix and the degree of freedom structures to be > initialized, and those two things are probably what's sucking up all > your memory. I see. > My goal is just to output the assembled system matrices to files and > > i don't have to solve them inside libmesh. > > The system matrix should have 551368 degrees of freedom, most of which > couple to 216 others. With 8 byte coefficients that's a hundred megs > of RAM and with sparsity pattern overhead it's probably two hundred > megs... but nine hundred MB seems excessive. Are you solving for more > than one scalar, using a system like a generalized EigenSystem that > builds more than one matrix, using complexvalued variables, or > anything else that might bump up the RAM requirements? It's for solving 2 system matrices in a realvalued generalized eigenvalue problem. I went back and tried (30x30x30, HEX27, 3rd order HERMITE). This time it didn't blew out the memory, but did take significantly longer time to assemble the matrices than (30x30x30, HEX27, 2nd order LANGRANGE). Maybe HERMITE is the problem? I'd appreciate it if you've got debugging tools that can give you a > memory breakdown by object type and you could give use such output. > It sounds like either we or PETSc might need to do a little more > optimization. I don't have any debugging tools and I have to admit that I'm a belowaverage C++ user. To work around your immediate problem, however: can you output the > element matrices instead of the system matrix, and assemble them > outside of libMesh? It sounds like you're using a structured mesh, > which can require much less overhead than the unstructured mesh class > in libMesh. How to assemble the element matrices outside of libmesh? Do you know any existing code/program that can do that? So, that would be, output each element matrix to a file and a problem should be able to read in all the element matrices from the file and assemble them to system matrices. I'm definitely interested if this is doable. Thanks, David 
From: Roy Stogner <roystgnr@ic...>  20060722 23:53:53

On Sat, 22 Jul 2006, David Xu wrote: > I was trying to assemble the stiffness and mass matrices on a dense mesh > with 531441 nodes (40x40x40, HEX27, 3rd, HERMITE) Can I suggest trying HEX8? The HERMITE elements are unique among our higher order elements in that all their degrees of freedom are topologically associated with mesh vertices, so unless you need quadratic mapping functions you don't need HEX27 elements. Those nodes aren't reponsible for most of your memory use (the system matrix is), but every megabyte helps. > and I ran into "Out of memory" problem at the line of: > > equation_systems.init(); > > > Here's the error message: > > [0]PETSC ERROR: PetscMallocAlign() line 62 in src/sys/memory/mal.c > [0]PETSC ERROR: Out of memory. This could be due to allocating > [0]PETSC ERROR: too large an object or bleeding by not properly > [0]PETSC ERROR: destroying unneeded objects. > [0]PETSC ERROR: Memory allocated 1380430780 Memory used by process 886095872 > [0]PETSC ERROR: Try running with malloc_dump or malloc_log for info. > [0]PETSC ERROR: Memory requested 907039528! > [0]PETSC ERROR: PetscTrMallocDefault() line 191 in src/sys/memory/mtr.c > [0]PETSC ERROR: MatSeqAIJSetPreallocation_SeqAIJ() line 2735 in > src/mat/impls/aij/seq/aij.c > [0]PETSC ERROR: MatCreateSeqAIJ() line 2621 in src/mat/impls/aij/seq/aij.c > [0]PETSC ERROR: User provided function() line 137 in > unknowndirectory/src/numerics/petsc_matrix.C > [unset]: aborting job: > application called MPI_Abort(comm=0x84000000, 1)  process 0 > > My question is: is it possible to assemble the matrices without having to > initilize the equation system? I'm afraid not. You could initialize a finite element object and evaluate the element matrices, but putting them into the system matrix requires that matrix and the degree of freedom structures to be initialized, and those two things are probably what's sucking up all your memory. > My goal is just to output the assembled system matrices to files and > i don't have to solve them inside libmesh. The system matrix should have 551368 degrees of freedom, most of which couple to 216 others. With 8 byte coefficients that's a hundred megs of RAM and with sparsity pattern overhead it's probably two hundred megs... but nine hundred MB seems excessive. Are you solving for more than one scalar, using a system like a generalized EigenSystem that builds more than one matrix, using complexvalued variables, or anything else that might bump up the RAM requirements? I'd appreciate it if you've got debugging tools that can give you a memory breakdown by object type and you could give use such output. It sounds like either we or PETSc might need to do a little more optimization. To work around your immediate problem, however: can you output the element matrices instead of the system matrix, and assemble them outside of libMesh? It sounds like you're using a structured mesh, which can require much less overhead than the unstructured mesh class in libMesh.  Roy 