the problem of radial

[Guangyuan Wei]

Dear developer,

I am considering phase transition in a two-dimensional disk using polar coordinates, with a constant pressure boundary condition. Since the pressure is a function of density, this is transformed into a constant density boundary condition. For the velocity boundary condition, I have applied a no-stress boundary condition. The equations are represented as follows:

******************************************************************************

k = (rho*d3.div(d3.grad(rho))+0.5*d3.grad(rho)@d3.grad(rho))*I-d3.grad(rho)*d3.grad(rho)

vpstress = -p*I + (1/Re)*(d3.grad(u) + d3.trans(d3.grad(u))-(2/3)*d3.div(u)*I)

kstress =  (1/We)*k

problem.add_equation("radial(vpstress)(r=1) = radial(-kstress)(r=1)")

problem.add_equation("(rho(r=1) = 0.4")

******************************************************************************

However, when I used the radial method to extract the normal stress, the following error occurred.

******************************************************************************

NotImplementedError: No subclasses of <class 'dedalus.core.operators.RadialComponent'> found for the supplied arguments: [Add(Convert(Mul(Mul(-1, <Field 140532601133184>), <Field 140532601127664>)), Mul(0.01, Add(Add(Grad(<Field 140532682543216>), TransposeComponents(Grad(<Field 140532682543216>))), Mul(-1, Mul(Mul(0.6666666666666666, Div(<Field 140532682543216>)), <Field 140532601127664>)))))], {'index': 0, 'out': None}

******************************************************************************

Additionally, I used the following method to obtain the total normal stress, but the numerical algorithm diverged. I would like to know if this method is correct. The code is shown as follows:

******************************************************************************

er = dist.VectorField(coords, bases=disk.radial_basis)

er['g'][0],  er['g'][1] = 0,  1

ephi = dist.VectorField(coords, bases=disk.radial_basis)

ephi['g'][0],  ephi['g'][1] = 1,  0
problem.add_equation("(er@vpstress)(r=1) = (-er@kstress)(r=1)")

******************************************************************************

How should I correct the code for this problem? I would be very grateful if you could help me.

Best wishes,

Guangyuan Wei

[Keaton Burns]

Hi Guangyuan,

The unit vector fields are not smooth fields on the disk, and this will likely cause some problems. You’re better off weighting them by r to get vector fields that are smooth, but still have the desired orientation at the boundary. Unfortunately we are still working on more flexible methods for applying component-wise boundary conditions in curvilinear domains, so taking inner products with vector fields is the best option at this time.

Best,

-Keaton

--
You received this message because you are subscribed to the Google Groups "Dedalus Users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to dedalus-user...@googlegroups.com.
To view this discussion visit https://groups.google.com/d/msgid/dedalus-users/3b554421-f257-4822-b386-dbede23df8c7n%40googlegroups.com.

[Guangyuan Wei]

Dear developer,

Thank you for your reply. I tried this a few days ago and the result was not satisfactory. So I changed the boundary conditions and wanted to ask you a few questions.

First I want to solve for the viscous stress term by representing the unit vector in polar coordinates. The equation is shown below:

I used the following representation, but I don't know if it's correct.

*******************************************************************************

I = dist.TensorField((coords, coords),bases=disk.radial_basis)

I['g'][0,0] = 1

I['g'][1,1] = 1

turu = (d3.grad(u) + d3.trans(d3.grad(u))-(2/3)*d3.div(u)*I)

*******************************************************************************

Second, I replace the stress-free boundary condition with the directional derivative along the normal direction, and the equation is as follows:

I used the following representation, but I don't know if it's correct.

*******************************************************************************

er = dist.VectorField(coords, bases=disk.radial_basis)

er['g'][0],  er['g'][1] = 0,  r

problem.add_equation("(er@grad(u))(r=1) = 0", condition='nphi!=0')

problem.add_equation("(er@grad(u))(r=1) = 0", condition='nphi==0')

problem.equations[-1]['LHS'].valid_modes[1] = True

*******************************************************************************

Is there a better expression for solving these two problems, making the numerical solution more stable? I would be grateful if you could reply.

Best wishes,

Guangyuan Wei