[fpena]

Dear Foamers,

I'm trying to use the phiformulation, im using a good solver named epotFoam created by Alessandro.

But im having some trouble when i'm trying to solve hartmann >0, i get the H=0 case any time, i made some modifications in the createFields to use the funkySetFields tool, and set a specific B-field

Info<< "Reading field B\n" << endl;
    volVectorField B0
    (
        IOobject
        (
            "B0",
            runTime.timeName(),
            mesh,
            IOobject::MUST_READ,
            IOobject::NO_WRITE
        ),
        mesh
    );

and i dropped the electricalProperties file in the constant directory.

Could anybody please correct me? Or there is something I am missing?

#include "fvCFD.H"
#include "pisoControl.H"

int main(int argc, char *argv[])
{
	#include "setRootCase.H"
	#include "createTime.H"
	#include "createMesh.H"

	pisoControl piso(mesh);

	#include "createFields.H"
	#include "initContinuityErrs.H"

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
	Info<< "\nStarting time loop\n" << endl;

	volVectorField lorentz = sigma *
	(
		((U ^ B0) ^ B0) -
		(fvc::grad(PotE) ^ B0)
	);

	while (runTime.loop())
	{
		Info<< "Time = " << runTime.timeName() << nl << endl;

		#include "CourantNo.H"

		fvVectorMatrix UEqn
		(
			fvm::ddt(U)
			+ fvm::div(phi, U)
			- fvm::laplacian(nu, U)
			==
			(1.0/rho) * lorentz
		);

		if (piso.momentumPredictor())
		{
			solve(UEqn == -fvc::grad(p));
		}

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
		while (piso.correct())
		{
			volScalarField rAU(1.0/UEqn.A());
			volVectorField HbyA(constrainHbyA(rAU*UEqn.H(), U, p));

			surfaceScalarField phiHbyA
			(
				"phiHbyA",
				fvc::flux(HbyA)
				+ fvc::interpolate(rAU)
				* fvc::ddtCorr(U, phi)
			);

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
			adjustPhi(phiHbyA, U, p);

			constrainPressure(p, U, phiHbyA, rAU);

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
			while (piso.correctNonOrthogonal())
			{

				fvScalarMatrix pEqn //Pressure equation coefficient matrix
				(
					fvm::laplacian(rAU, p)
					==
					fvc::div(phiHbyA)
				);

				pEqn.setReference(pRefCell, pRefValue);

				pEqn.solve(mesh.solver(p.select(piso.finalInnerIter())));

				if (piso.finalNonOrthogonalIter())
				{
					phi = phiHbyA - pEqn.flux();
				}
			}

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
			#include "continuityErrs.H"

			U = HbyA - rAU*fvc::grad(p);
			U.correctBoundaryConditions();
		}

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
		surfaceScalarField psiub = fvc::interpolate(U ^ B0) & mesh.Sf();

		fvScalarMatrix PotEEqn
		(
			fvm::laplacian(PotE)
			==
			fvc::div(psiub)
		);

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
		PotEEqn.setReference(PotERefCell, PotERefValue);

		PotEEqn.solve();

		surfaceScalarField jn = psiub - (fvc::snGrad(PotE) * mesh.magSf());

		surfaceVectorField jnv = jn * mesh.Cf();

		volVectorField jfinal = fvc::surfaceIntegrate(jnv) - (fvc::surfaceIntegrate(jn) * mesh.C());

		jfinal.correctBoundaryConditions();

		lorentz = sigma* (jfinal ^ B0);

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
		runTime.write();

		Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
		<< " ClockTime = " << runTime.elapsedClockTime() << " s"
		<< nl << endl;
		//end of time loop
	}

	Info<< "End\n" << endl;
	return 0;
}

(i will continue traying by myselft oc); i attach the solver epotFoam, and the case how is set.

Appreciate your time and effort in advance!

Best regards.