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}
(62) 
_{Sep}

_{Oct}

_{Nov}

_{Dec}

S  M  T  W  T  F  S 







1

2

3

4
(1) 
5
(1) 
6

7
(8) 
8

9

10
(5) 
11
(6) 
12
(2) 
13
(1) 
14
(1) 
15
(1) 
16

17

18

19

20

21

22

23

24

25
(1) 
26

27

28

29

30

31






From: John Peterson <peterson@cf...>  20050111 16:16:29

Mark Blome writes: > > Hi John and Michael, > > thanks for your comments on this problem! I think I finally found the mistake. > As usually in a place where I never expected it: The right hand side terms > where correct, but for the calculation of the gradient of the potential Un I > used something like > > double x= qPoints[qp](0)source.pos(0) ; > double y= qPoints[qp](1)source.pos(1) ; > double z= qPoints[qp](2)source.pos(2) ; > > double r=sqrt( pow(x,2) + pow(y,2) + pow(z,2) ); > double cpart=  source.j0/(2*libMesh::pi * sigma0); > RealGradient gradVn; > gradVn(0)=cpart*x/pow(r,3); gradVn(1)=cpart*y/pow(r,3); > gradVn(2)=cpart*z/pow(r,3); > > in the quadrature loop. > (source.pos is the coordinate of the current source and the qpoints are from > the quad rule: const std::vector<Point> qPoints=qrule.get_points();) > I was simply not aware that "get_points()" returns local coordinates on a > reference element and not global coodinates. If I use eg "elem>point(0)" > instead of qPoints[qp] the results are fine! Is there an easy way to convert > the local reference element coordinates to global coordinates? Yes, use get_xyz(). John 
From: Mark Blome <mblome@au...>  20050111 16:06:57

Hi John and Michael, thanks for your comments on this problem! I think I finally found the mistake. As usually in a place where I never expected it: The right hand side terms where correct, but for the calculation of the gradient of the potential Un I used something like double x= qPoints[qp](0)source.pos(0) ; double y= qPoints[qp](1)source.pos(1) ; double z= qPoints[qp](2)source.pos(2) ; double r=sqrt( pow(x,2) + pow(y,2) + pow(z,2) ); double cpart=  source.j0/(2*libMesh::pi * sigma0); RealGradient gradVn; gradVn(0)=cpart*x/pow(r,3); gradVn(1)=cpart*y/pow(r,3); gradVn(2)=cpart*z/pow(r,3); in the quadrature loop. (source.pos is the coordinate of the current source and the qpoints are from the quad rule: const std::vector<Point> qPoints=qrule.get_points();) I was simply not aware that "get_points()" returns local coordinates on a reference element and not global coodinates. If I use eg "elem>point(0)" instead of qPoints[qp] the results are fine! Is there an easy way to convert the local reference element coordinates to global coordinates? Mark Am Dienstag 11 Januar 2005 01:29 schrieb John Peterson: > Hi Mark, > > Can you post a compilable version of your code? My only other > guess is you could have a sign error since it's going to the right > hand side. > > J > > Mark Blome writes: > > Hi Michael, > > > > I am using libmesh to solve for the electric potential of a current > > source within a 3d cubic domain. The calculations are used for a > > geoelectrical forward solver. On the top side of the domain, I apply a > > Neuman boundary condition (du/dn=0) and on the other sides I use > > Dirichlet boundary conditions (u=0). > > The second term > > div( (sigmasigma0) grad Un) > > is added to remove the singularity arising due to the point current > > source. When I discretize the second term with > > > > Fe(i) +=  (sigmasigma0) * JxW[qp] * (gradVn[qp] * dphi[i][qp]) > > (of course, as John noticed, I forgot to mention the surface integral > > term, but in my case it is zero anyway cause sigma=sigma0 at the outer > > boundary. ) > > > > I get wrong results: The resulting potential U is the potential that > > does not contain the singular potential Un=1/(2*pi) *1/r anymore. > > For the case of a homogenous cube with a high conducting cube centered > > inside, the introduction of the discretization above results in a > > strange dipole field within the high conducting cube, which physically > > doesnt make any sense. > > Thanks for the help so far, > > > > Mark > > > > Am Montag 10 Januar 2005 19:28 schrieb Michael Schindler: > > > Hi Mark, > > > > > > On 10.01.05, Mark Blome wrote: > > > > Hi everybody, > > > > > > > > I am trying to descretize a laplace equation which contains a second > > > > term. In my opinion the problem should be very easy, but I become > > > > more and more confused the more I think about it. The laplace > > > > equation reads: div(sigma grad U) + div( (sigmasigma0) grad Un) =0 > > > > where Un is a potential for which the equation is known ( Un=j0/(2 * > > > > PI * sigma0) * 1/r), hence the term has to go to the right hand > > > > side. > > > > > > > > For the term div(sigma grad U)=0 I simply get matrix terms: > > > > Ke(i,j) +=  sigma * JxW[qp]*(dphi[i][qp]*dphi[j][qp]); > > > > But what to add for div( (sigmasigma0) grad Un)? > > > > > > The FEM method is to find a weak solution with a testfunction phi > > > > > > int div(sigma grad U) phi = int div( (sigmasigma0) grad Un) phi > > > > > > Depending on sigma you can use the Galerkinansatz where phi is the > > > same function you approximate U with. If sigma is constant, omitting > > > the surface terms. > > > > > >  sigma int grad(U)*grad(phi) =  (sigmasigma0) int > > > grad(Un)*grad(phi) > > > > > > This is what you tried first. > > > > > > > I tried: > > > > Fe(i) +=  (sigmasigma0) * JxW[qp] * (gradVn[qp] * > > > > dphi[i][qp]) ; where gradVn[qp] simply is grad(Un) analytically > > > > calculated at the quadrature points. But it doesnt work like this. > > > > > > What actually does not work? What boundary conditions do you use? > > > > > > Michael. > > > >  > > The SF.Net email is sponsored by: Beat the postholiday blues > > Get a FREE limited edition SourceForge.net tshirt from ThinkGeek. > > It's fun and FREE  well, almost....http://www.thinkgeek.com/sfshirt > > _______________________________________________ > > Libmeshusers mailing list > > Libmeshusers@... > > https://lists.sourceforge.net/lists/listinfo/libmeshusers 
From: John Peterson <peterson@cf...>  20050111 14:40:14

Martin L=FCthi writes: > Hi >=20 > Maybe a stupid question: Why is the following construct used in all > the examples to determine that an element lives on the boundary? >=20 > for (unsigned int s=3D0; s<elem>n_sides(); s++) > if (elem>neighbor(s) =3D=3D NULL) >=20 > Would'nt this be faster / more obvious if one would use >=20 > if (elem>on_boundary()) >=20 > and then use the boundary side iterators (Elem::boundary_sides_begin= )? Heh, yes probably. The boundary side iterators were only recently added, the rest of the code hasn't been updated yet. John 
From: KIRK, BENJAMIN (JSCEG) (NASA) <benjamin.kirk1@na...>  20050111 14:30:37

Absolutely. This is just a case of the examples predating this method = and the boundary side iterators. They should be updated. Along those lines, I have a new approach for building sides that is = less expensive in terms of memory allocation and data copying which will be included in libMesh0.5.0 realsoonnow. Ben Original Message From: libmeshusersadmin@... [mailto:libmeshusersadmin@...] On Behalf Of Martin = L=FCthi Sent: Monday, January 10, 2005 7:28 PM To: libmeshusers@... Subject: [Libmeshusers] Testing if element is on boundary Hi Maybe a stupid question: Why is the following construct used in all the examples to determine that an element lives on the boundary? for (unsigned int s=3D0; s<elem>n_sides(); s++) if (elem>neighbor(s) =3D=3D NULL) Would'nt this be faster / more obvious if one would use if (elem>on_boundary()) and then use the boundary side iterators (Elem::boundary_sides_begin)? Thanks for clarification Martin  The SF.Net email is sponsored by: Beat the postholiday blues Get a = FREE limited edition SourceForge.net tshirt from ThinkGeek. It's fun and = FREE  well, almost....http://www.thinkgeek.com/sfshirt _______________________________________________ Libmeshusers mailing list Libmeshusers@... https://lists.sourceforge.net/lists/listinfo/libmeshusers 
From: <luthi@gi...>  20050111 01:22:32

Hi Maybe a stupid question: Why is the following construct used in all the examples to determine that an element lives on the boundary? for (unsigned int s=0; s<elem>n_sides(); s++) if (elem>neighbor(s) == NULL) Would'nt this be faster / more obvious if one would use if (elem>on_boundary()) and then use the boundary side iterators (Elem::boundary_sides_begin)? Thanks for clarification Martin 
From: John Peterson <peterson@cf...>  20050111 00:30:07

Hi Mark, Can you post a compilable version of your code? My only other guess is you could have a sign error since it's going to the right hand side. J Mark Blome writes: > > Hi Michael, > > I am using libmesh to solve for the electric potential of a current source > within a 3d cubic domain. The calculations are used for a geoelectrical > forward solver. On the top side of the domain, I apply a Neuman boundary > condition (du/dn=0) and on the other sides I use Dirichlet boundary > conditions (u=0). > The second term > div( (sigmasigma0) grad Un) > is added to remove the singularity arising due to the point current source. > When I discretize the second term with > > Fe(i) +=  (sigmasigma0) * JxW[qp] * (gradVn[qp] * dphi[i][qp]) > (of course, as John noticed, I forgot to mention the surface integral term, > but in my case it is zero anyway cause sigma=sigma0 at the outer boundary. ) > > I get wrong results: The resulting potential U is the potential that does not > contain the singular potential Un=1/(2*pi) *1/r anymore. > For the case of a homogenous cube with a high conducting cube centered inside, > the introduction of the discretization above results in a strange dipole > field within the high conducting cube, which physically doesnt make any > sense. > Thanks for the help so far, > > Mark > > > > Am Montag 10 Januar 2005 19:28 schrieb Michael Schindler: > > Hi Mark, > > > > On 10.01.05, Mark Blome wrote: > > > Hi everybody, > > > > > > I am trying to descretize a laplace equation which contains a second > > > term. In my opinion the problem should be very easy, but I become more > > > and more confused the more I think about it. The laplace equation reads: > > > div(sigma grad U) + div( (sigmasigma0) grad Un) =0 > > > where Un is a potential for which the equation is known ( Un=j0/(2 * PI * > > > sigma0) * 1/r), hence the term has to go to the right hand side. > > > > > > For the term div(sigma grad U)=0 I simply get matrix terms: > > > Ke(i,j) +=  sigma * JxW[qp]*(dphi[i][qp]*dphi[j][qp]); > > > But what to add for div( (sigmasigma0) grad Un)? > > > > The FEM method is to find a weak solution with a testfunction phi > > > > int div(sigma grad U) phi = int div( (sigmasigma0) grad Un) phi > > > > Depending on sigma you can use the Galerkinansatz where phi is the > > same function you approximate U with. If sigma is constant, omitting > > the surface terms. > > > >  sigma int grad(U)*grad(phi) =  (sigmasigma0) int grad(Un)*grad(phi) > > > > This is what you tried first. > > > > > I tried: > > > Fe(i) +=  (sigmasigma0) * JxW[qp] * (gradVn[qp] * dphi[i][qp]) > > > ; where gradVn[qp] simply is grad(Un) analytically calculated at the > > > quadrature points. But it doesnt work like this. > > > > What actually does not work? What boundary conditions do you use? > > > > Michael. > > >  > The SF.Net email is sponsored by: Beat the postholiday blues > Get a FREE limited edition SourceForge.net tshirt from ThinkGeek. > It's fun and FREE  well, almost....http://www.thinkgeek.com/sfshirt > _______________________________________________ > Libmeshusers mailing list > Libmeshusers@... > https://lists.sourceforge.net/lists/listinfo/libmeshusers  John W. Peterson The University of Texas at Austin Web: http://www.cfdlab.ae.utexas.edu Phone: 5126196534 (H)  5124714069 (W) email: peterson@... 