Regarding magnetic fluid behaviour under non-uniform magnetic field

[Shyam Sunder Yadav]

Dear Prof. Popinet and Prof. Lopez

I am trying to simulate flow behavior of a magnetic fluid under a non-uniform magnetic field assuming linear magnetization behavior.

The equations governing the magnetic field related phenomenon are as follows:

I have written a very simple code to implement the above equations.

I am applying equal and opposite magnetic potential to the top and bottom boundaries of the domain.

The distribution of magnetic potential, magnetic field components is correct and as expected.

The issue I am facing is that the resulting velocity field is not coming symmetric at the top and bottom boundaries.

Please check the attached movie and screenshot.

I have attached the code also. Please let me know what mistake I am making.

Particularly, please check whether the implementation of acceleration event in the file magnet.h is correct or not. 

Thank you. Regards

--

Dr. Shyam Sunder Yadav
Associate Professor
Mechanical Engineering
BITS Pilani
09902346342
http://www.bits-pilani.ac.in/pilani/ssyadav/Profile

The information contained in this electronic communication is intended solely for the individual(s) or entity to which it is addressed. It may contain proprietary, confidential and/or legally privileged information. Any review, retransmission, dissemination, printing, copying or other use of, or taking any action in reliance on the contents of this information by person(s) or entities other than the intended recipient is strictly prohibited and may be unlawful. If you have received this communication in error, please notify us by responding to this email or telephone and immediately and permanently delete all copies of this message and any attachments from your system(s). The contents of this message do not necessarily represent the views or policies of BITS Pilani.

[Shyam Sunder Yadav]

Dear Sir

I tried the Maxwell stress tensor based formulation as well adapted from ehd/stress.h 

Still the velocity field is unsymmetric about the horizontal axis...

I think physically it should be symmetric for linear magnetization and at initial times...

Please look into this.

Thanks and Regards

[Robert A]

Dear  Shyam,

Adding the following boundary conditions for the face fluid velocity:

uf.n[left]   = dirichlet(0.);
uf.n[right]  = dirichlet(0.); //neumann(0.);
uf.n[bottom] = dirichlet(0.);
uf.n[top]    = dirichlet(0.);

uf.t[left]   = dirichlet(0.);
uf.t[right]  = dirichlet(0.); //neumann(0.);
uf.t[bottom] = dirichlet(0.);
uf.t[top]    = dirichlet(0.);

does make the solution perfectly symmetric about the x-axis (as you rightfully stated).

I am myself surprised by the fact that without these bc's the solution is not symmetric even for a perfectly symmetric acceleration/force term a[]. 

All the best,

Robert 

[Edoardo Cipriano]

Dear Shyam and Robert,

I am not 100% sure, but I think that the boundary condition on the normal component of uf should not use the function `dirichlet()`. Instead, you directly set the value, for example:

uf.n[top] = 0.;

You can verify how this is done in the available Basilisk codes using the commands:

grep -r 'uf.n\[' ~/basilisk/src/

grep -r 'uf.t\[' ~/basilisk/src/

Best,

Edoardo

[Stephane Popinet]

> I think that the boundary condition on the normal component of uf
> should not use the function `dirichlet()`. Instead, you directly set
> the value, for example:
>
> uf.n[top] = 0.;

This is correct.

Stephane

[Robert A]

Thanks, Stephane. Good to know.

Many cases in the sandbox with dirichet(0.) will then have to be revised (see attached image for example)

[Robert A]

Dear Edoardo,

Setting

uf.n[left]   = 0.;
uf.n[right]  = 0.;
uf.n[bottom] = 0.;
uf.n[top]    = 0.;

worked too. (AFAIK this sets boundary values directly in the ghost cells as opposed to those on the boundary faces)

The question is: why do we have to set boundary conditions on uf at all if vast majority of Basilisk setups do not involve bc's on uf. (There are some exceptions, see attached)

[Stephane Popinet]

Correction:

uf.n[left] = 0.;

is correct, but

uf.n[left] = dirichlet(0.);

is also correct (and will do exactly the same thing).

Stephane

[Shyam Sunder Yadav]

Dear Robert, Edoardo, Prof. Popinet 

Thank you very much for the information!

The velocity field is indeed symmetric with the face based boundary conditions for velocity.

This is the case for both Kelvin body force as well as divergence of Maxwell stress tensor based formulations.

But the flow field pattern (vortex advection speed, pattern etc.) is different for Kelvin force and Maxwell stress tensor based formulations.

Things also depend on whether we are storing magnetic field components on faces or cell centers.

 I will discuss these issues in the coming September BMM.

Thanks again dear all for the suggestions.

Regards

--
You received this message because you are subscribed to the Google Groups "basilisk-fr" group.
To unsubscribe from this group and stop receiving emails from it, send an email to basilisk-fr...@googlegroups.com.
To view this discussion visit https://groups.google.com/d/msgid/basilisk-fr/c2d2ebd9-3268-414e-a86d-f0580b67cc8b%40basilisk.fr.

[jose.lopez...@gmail.com]

Hi,

uf is a face magnitud constructed by interpolating  the velocity from cell centered values(i.e uf.x[] = (v.x[-1] + v.x[])/2) and fundamentally serves to compute convective fluxes into the cell. Doing uf.n[left]=0 you are forcing the flux value at the border of the computational domain directly obviating the interpolation (consequently the ghost cell value is not used for the computation of fluxes at that border)

Cheers

pepe

[Shyam Sunder Yadav]

Dear Prof. Lopez

Thank you very much for the insight...

Regards

--

You received this message because you are subscribed to the Google Groups "basilisk-fr" group.
To unsubscribe from this group and stop receiving emails from it, send an email to basilisk-fr...@googlegroups.com.

To view this discussion visit https://groups.google.com/d/msgid/basilisk-fr/da1d0e36-e281-4756-abc5-d53358ac6e74n%40googlegroups.com.

[Robert A]

Glad that it worked out for you, Shyam.

If you look at the force components a.x and a.y near x=-0.02 && x<0.02 close the top and bottom boundaries, they jump strongly from one cell to another.

You can try to make phi vary from 0 to 1000 and back to 0 smoothly at the boundary instead of the step-like boundary condition. 

Then differences between storing magnetic field components on faces or cell centers may become less prominent. 

You can also try to activate cell refinement on u.x and u.y.

All the best,

Robert

[Robert A]

Hi Pepe,

Thanks for your clarification!

If you set u.n[left] = dirichlet(0.) this should mean uf.n[left] = 0 automatically?

If say, u.x = 5 in a cell attached to a left boundary cell, they the ghost cell value will be -5, so that (5+(-5))/2 = 0 on the face?

All the best,

Robert

[Robert A]

Glad that it worked for you, Shyam.

It's still not very clear to me when to apply and not apply boundary conditions for uf.

You may want to look to the force term components a.x and a.y around x=0.02 and x = -0.02 at the top and bottom boundaries. They vary a lot from one cell to another. 

If you can formulate a smoother boundary condition for phi, then storing magnetic field components (needed only to compute a[]) on faces or cell centers may not matter much. 

Of you can try to refine the grid around x=0.02 and x = -0.02.

All the best, 

Robert

On Thursday, July 31, 2025 at 12:54:43 PM UTC-5 jose.lopez...@gmail.com wrote:

[Victor Boniou]

Hi everyone, 

This might be off-topic, but I think that you have a typo in your acceleration event:

foreach_face() <------- Should be foreach_face(x) otherwise you update twice av.y[] with possibly wrong values! { // face values of magnetic field components double H_x = H.x[]; double H_y = 0.25*(H.y[]+H.y[0,1]+H.y[-1,1]+H.y[-1,0]); // f = mu0 * Chi * H \cdot \nabla H av.x[] += (alpha.x[]/fm.x[]) * Mu0 * Chi * ( H_x*(H.x[1,0]-H.x[-1,0]) + H_y*(H.x[0,1]-H.x[0,-1]) )/(2.0*Delta); } /* foreach_face(y) { double H_y = H.y[]; double H_x = 0.25*(H.x[]+H.x[0,-1]+H.x[1,-1]+H.x[1,0]); av.y[] += (alpha.y[]/fm.y[]) * Mu0 * Chi * ( H_x*(H.y[1,0]-H.y[-1,0]) + H_y*(H.y[0,1]-H.y[0,-1]) )/(2.0*Delta); }

Hope this can help you,

Best,

Victor

[Shyam Sunder Yadav]

Dear Victor

Thanks for the email...

The  foreach_face(y){}  loop is commented out...

We correctly have the   foreach_face(){}  loop.

Please check again.

Thanks and Regards

--

You received this message because you are subscribed to the Google Groups "basilisk-fr" group.
To unsubscribe from this group and stop receiving emails from it, send an email to basilisk-fr...@googlegroups.com.

To view this discussion visit https://groups.google.com/d/msgid/basilisk-fr/9ce966af-2e2a-4684-a344-5ffe989d8207n%40googlegroups.com.