Mesh in polar coordinates

General Discussion
1

Suppose we have a domain in polar coordinates with boundary r=f(theta). 0<=theta<=2*pi. How to build a mesh in FreeFem++. Can I use movemesh in polar to polar transformation from (r , theta) coordinates to (r’ , theta’) coordinates. This FreeFem++ is far better than MatLab for FEMians. Any suggestion regarding this will be a great help.
Thank you!

2

Just do like a circle mesh:

func real f(real t){ return 1+0.1*cos(t);}
border bf(t=0,2*pi){ x=f(t)*cos(t);y=f(t)*sin(t}:
mesh Th=buildmesh(bf(30));
plot(Th,wait=1);
3

What about transformation of mesh from polar (r , theta) to polar (r’ , theta’). How to use movemesh ?

4

pass trough x,y coordinate

5

Monsieur,

Can we handle the polar coordinates using the functions such as written below?

//FE Velocity space
fespace Vh(Th,P2); // P2 elements

Vh u,v, uh, vh;

real eps = 1.0e-10;

func f1 = x;
func f2 = cos(y);
func f3 = sin(y);

// variational formulation
problem lubrication([u,v,p], [uh,vh,ph]) =
   int2d(Th) ( f1^2*nu*dx(u)*dx(uh) + f2^2*nu*dy(u)*dy(uh) )
 + int2d(Th) (f1^2*nu*dx(v)*dx(vh) + f2^2*nu*dy(v)*dy(vh) )
 - int2d(Th) (f1*p*dx(uh))
 - int2d(Th) (f2*p*dy(vh))
 - int1d(Th, Sh) (f1*nu*dx(u)*N.x*uh+f2*nu*dy(u)*N.y*uh)
 - int1d(Th, Ho) (f1*nu*dx(v)*N.x*vh+f2*nu*dy(v)*N.y*vh)
 + int2d(Th) (f1*dx(u)*ph + f2*dy(v)*ph)
 + int2d(Th) (eps*p*ph)
  + on (1, u=0, v=V0)
 + on (2, u=0, v=1.2);

lubrication;

I am new to coding in FreeFEM++ and borrowed the idea of generating a 3D cylinder from cube. I don’t know the correctness. But, I got a better result from this code, so a bit surprised.

Thanks and best regards,

6

Sorry, Th is a square or a rectangle, I do not understand what you do.
I need the full script to be sure.

7

Dear Sir,

Thank you for your response. Please see the code for the journal bearing and shaft arrangement here -

load "msh3"
load "shell"
load "medit"

// Body geometry

// Mesh quality
int nu1 = 50;  // Mesh quality for the bearing shell
int nu2 = 200; // Mesh quality for the obstacle
int nu3 = 500;  // Mesh quality for the shaft
int nu4 = 500; // Mesh quality for the lubricant
real ldratio=2; // Bearing L/diameter ratio

// Dimensions of the shaft, bearing and bearing housing
real Ri=0.350; // Shell bearing inner radius
real Ro=0.380; // Shell bearing outer radius
real clr=0.010 + 3e-4; // Shaft - Bearing clearance
real eccx=0.0, eccy=-0.0100; // Eccentricity of the shaft in bearing housing
real Rs=Ri-clr; // Shaft radius
real wear = 0.0001; // bearing uniform wear
real thck = Ro-Ri; // Bearing thickness

//Propeller shaft
border Sh(t=0,2*pi){x=Rs*cos(t)+eccx; y=Rs*sin(t)+eccy; label=0;}
mesh Shft = buildmesh(Sh(nu1));

//Top Bearing Shell
border bt1(t=pi,0){x=Ri*cos(t);y=Ri*sin(t);label=1;}
border bt2(t=0,thck){x=Ri+t;y=0.;label=2;}
border bt3(t=0,pi){x=Ro*cos(t);y=Ro*sin(t);label=3;}
border bt4(t=0,thck){x=-Ro+t;y=0;label=4;}
mesh Thbt=buildmesh(bt1(nu1) +bt2(nu1) +bt3(nu1) +bt4(nu1));
//plot(Th1t, wait=true);

//Bottom Bearing Shell
border bb1(t=0,pi){x=Ri*cos(t);y=Ri*sin(-t);label=11;}
border bb2(t=0,thck){x=-Ri-t;y=0.;label=12;}
border bb3(t=pi,0){x=Ro*cos(t);y=Ro*sin(-t);label=13;}
border bb4(t=0,thck){x=Ro-t;y=0;label=14;}
mesh Thbb=buildmesh(bb1(nu1) +bb2(nu1) +bb3(nu1) +bb4(nu1));
//plot(Th1b, wait=1);

//Lubricant Domain
mesh Thlub = buildmesh(bb1(-nu4)+bt1(-nu4)+Sh(-nu3));
//Thlub = adaptmesh(Thlub,eccx*eccy);
//plot(Thlub, wait=1);


// Bearing Housing mesh
real R1=Ri+0.07;
real e2=R1-Ro;
border b21(t=0,2*pi){x=R1*cos(t);y=R1*sin(t);label=101;}
border b22(t=0,e2){x=-R1+t;y=0.;label=102;}
border b23(t=2*pi,0){x=Ro*cos(t);y=Ro*sin(t);label=103;}
border b24(t=0,e2){x=Ro+t;y=0;label=104;}
mesh Th1=buildmesh(b21(nu2)+b22(nu2/10)+b23(nu2)+b24(nu2/10));

mesh Th2 = Th1+Shft+Thbt+Thbb+Thlub;

mesh3 Th = buildlayers(Th2, 20, zbound=[0,ldratio*Ri]);

plot(Th2, wait=1);

medit("Stern Bearing Assy", Th);

/////////////////////////////////////////////////////////////////////////////////////////////////////
/*   Polar coordinates guidance                                                                  */
/*   x = r*cos(theta)                                                                                       */
/*   y = r*sin(theta)                                                                                        */
/*   (x,y) ----------> (r,theta)                                                                              */
/*                                                                                                                      */
/*   dr      -------->  cos(theta) * dx + sin(theta) * dy                                    */
/*   d(theta)  -------->  1/r * (-dx*sin(theta) + dy*cos(theta) )                    */
/*   dr * d(theta) ---->  1/r * dx*dy                                                                */
/*                                                                                                                      */
/*  cos(theta) = x/sqrt(x^2 + y^2)                                                                */
/*  sin(theta)  = y/sqrt(x^2 + y^2)                                                                */
/*                                                                                                                      */
/////////////////////////////////////////////////////////////////////////////////////////////////////

There are four separate regions for simulating the running of a shaft in a bearing housing with thin shell bearings.

  1. Shaft
  2. Lubricant film in the clearance region
  3. Top and bottom shells (bearings)
  4. Bearing housing

The variational formulation is for the lubricant film in the clearance region to determine the radial, whirl velocities and hydrodynamic pressure buildup as the shaft rotates at an RPM.

I want to use the Navier Stokes Equations in polar coordinates for the regions with the polar coordinates functions indicated in the code and solve. Is this the right approach?