[Top][All Lists]

[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

[ff3d-users] Problem with elliptic hollow cylinder section

From: Tobias Nähring
Subject: [ff3d-users] Problem with elliptic hollow cylinder section
Date: Thu, 22 Apr 2004 21:44:13 +0200
User-agent: Mozilla/5.0 (X11; U; Linux i686; en-US; rv:1.6) Gecko/20040401 Debian/1.6-4

A short summary of the problem with the corresponding pictures and
files is available at: http://www.tn-home.de/Soft/Fem3d/index.html.

I want to solve the Laplacian equation in the domain

 G := {(x,y,z)|   1 < (x/rx)^2 + (y/ry)^2 < a^2;   0 < z < L }

which resembles a section of a hollow cylinder with elliptical
cross-section. The constants are: rx=1, ry=0.5, a=2, L=5.

The boundary condition is

 u(x,y,z) = (|x,y|-1)*(a-|x,y|)*(L-z)/L


 |x,y|:=sqrt( (x/rx)^2 + (y/ry)^2 ).

That means effectively that u(x,y,z)=0 for z=L, and for (x,y,z) on the
cylinder mantles. Only at z=0 one has some nice non-zero boundary

Because of the symmetry of the boundary condition the
computational domain can be confined to the quarter

Q:={(x,y,z) in G| x>0, y>0}

of G with natural boundary conditions for u on the boundaries at x=0
and y=0.

I tried to solve this problem with FreeFEM3D without much luck.

The results of my attempts are shown at

ff3d accepts the input file and the corresponding pov file (included below)
and the solver options are set to the save values

cg(epsilon=1e-8) and cg(maxiter=2000)

but the output is rather erroneous. I have played
already with other solver settings, such as bicg(...) and
krylov(type=bicg) without any luck. The grid size is as large as
possible on my computer with 512MB RAM. (If I enlarge the grid size
then ff3d freezes during "Marching cube ...".

With the current grid size I get the error message:

 warning : meshing object plane {
 (0, 1, 0),2
  a problem occured ...
 Check that:
 - the object is not to small compare to the background mesh
 - the object intersects the background mesh

I don't know how to handle that.

The ff3d computed solution of the problem are graphically represented
in the figures:


Is there any chance to get a better solution?

One more remark: There would be a way to specify the mantle surface of
a cylinder section: the `open' directive of the povray object
cylinder (as described in the povray manual).

The input file for ff3d:
// Model of a cable with ellipsoidal cross-section.

vector Grid=(70,60,30);
double rx = 1.0; // radius of the inner conductor in x-direction
double ry = 0.5; // radius of the inner conductor in y-direction
double a=2.0; // The cross-section boundary of the inner conductor
         // scaled by a gives the cross-section boundary of the
         // shielding.
double CableW = a*rx;
double CableH = a*ry;

double CableL = 5.0; // length of the cable
double Border = 0.1; // Border added around the modelled piece of the
            // cable to build a fictitious domain

scene myScene = pov("cable-ellipse-annulus-num.pov");

mesh myMesh = structured( Grid, (-Border,-Border,-2.0*Border), (a*rx+Border,a*ry+Border,CableL+2.0*Border)); // Comment: in z-direcion we provided a bit more space for the fictitious domain (just guessing).

domain myDomain = domain(myScene, inside(<1,0,0>));
mesh myCylMesh = surface(myDomain,myMesh);

function norm2d=sqrt((x/rx)^2+(y/ry)^2);

// Realizing the boundary conditions with method-type penalty doesnot
// work (no convergence).
solve(u) in myDomain by myMesh
 pde(u) div(grad(u))=0;
 u=(norm2d-1)*(a-norm2d)*(CableL-z)/CableL on myCylMesh;

// Results output for postprocessing via gnuplot:

double myx; double myy; double myz;
double mydx = CableW/x(Grid);
double mydy = CableH/y(Grid);
double mydz = CableL/z(Grid);
double myxmax=CableW+0.5*mydx;
double myymax=CableH+0.5*mydy;
double myzmax=CableL+0.5*mydz;

for( myz=0.0; myz <= myzmax; myz = myz + mydz ) {
 for( myx=0.0; myx <= myxmax; myx = myx + mydx ) {
   for( myy=0.0; myy <= myymax; myy = myy + mydy ) {
cout << "@" << myx << ", " << myy << ", " << myz << ", " << u(myx,myy,myz) << "\n";
   cout << "@\n";
 cout << "@\n";

The pov file (named cable-ellipse-annulus-num.pov):
// Model of a cable with ellipsoidal cross-section.
// rx=1.0 radius of the inner conductor in x-direction
// ry=0.5 radius of the inner conductor in y-direction
// a=2.0 the inner conductor scaled by a gives the shielding
// CableL=5 length of the cable.
// Since we impose symmetric boundary conditions we need only
// to model a quarter of the cable. (0 < x < a*rx, 0 < y < a*ry)
intersection {
 difference {
   cylinder {
     2.0 // a==2.0
     <0,1,0>, 0
     <1,0,0>, 0
 scale <1.0,0.5,5.0> // <rx,ry,CableL>
 pigment { color rgb <1,0,0> }

reply via email to

[Prev in Thread] Current Thread [Next in Thread]