[Lang]

Hello everyone,

I am a newbie to OpenFoam programming. Currently, there is a challenge I'm facing in building a two-material continuous casting simulation (picture on the first page of this PDF for reference https://kupfer.de/fileadmin/user_upl...FA_Mittler.pdf)

The simulation involves four phases: two liquid phases (melts) and two solid phases. I initially attempted to use the "icoReactingMultiphaseInterFoam" solver, but I encountered an issue. The solver doesn't support moving the solid phase, which is crucial for pulling the strand. By browsing through the source code I found that the available phase models are:
- pureMovingPhaseModel
- pureStaticPhaseModel
- pureStaticSolidPhaseModel
- multiComponentMovingPhaseModel.
All of them do not meet my requirements. What I would need is a "pureMovingSolidPhaseModel".
I considered blending the available phase models to achieve a "pureMovingSolidPhaseModel," but I realized that it led me into a complex web of classes that I struggle to navigate in.

During my research, I came across another OpenFOAM solver called directChillFoam (https://github.com/blebon/directChillFoam) that specifically simulates a DC continuous casting process. In this solver, the strand pulling is simulated by introducing the velocity field of the solid (Us), without the need for a dynamic/moving mesh. This is exactly what I want to achieve in my simulation.

My idea was to modify the icoReactingMultiphaseInterFoam solver based on the approach used in directChillFoam. However, I'm unsure if this is feasible within a one-month timeframe, considering that directChillFoam is for OpenFOAM 9, which adds an additional challenge as I want to implement this in OF v2212.

I would greatly appreciate any tips, suggestions, or alternative approaches that you may have. Is it possible to create a new phase model from scratch to solve this problem? Or is there a way to achieve my goal using the icoReactingMultiphaseInterFoam solver that I haven't discovered yet?

Thank you in advance for your help and for taking the time to read my post!

Best regards,
Matthias

[Lang]

Greetings everyone,
I'd like to share the solution I've discovered for the issue I was facing. It turned out that the phaseModels were misleading me, leading me in the wrong direction. The actual problem stemmed from the implementation of the VollerParkash interface porosity model.
This model serves the purpose of simulating the "mushy zone" in a solidifying liquid. To achieve this, it introduces a source term on the diagonal of matrix A (referring to the linear system of equations Ax=b). This source term increases as the alpha.solid parameter grows, and it becomes extremely large (practically approaching infinity) when alpha.solid reaches 1. As a result, it enforces a velocity of zero in that specific cell.
The specific part of the source code responsible for this can be located in the file named multiphaseInterSystem.C. The method of interest is void Foam::multiphaseInterSystem::addInterfacePorosity( fvVectorMatrix& UEqn). In order to address the issue, I made modifications to the source code. Specifically, I adjusted the code so that the value added to the matrix diagonal is also added to the vector b (which is part of the Ax=b equation). This addition is multiplied by the desired solid velocity, which is read from a dictionary file.

void Foam::multiphaseInterSystem::addInterfacePorosity(fvVectorMatrix& UEqn)
{
    const scalarField& Vc = mesh_.V();
    scalarField& Udiag = UEqn.diag();
    vectorField& Usource = UEqn.source();

    forAllConstIters(phaseModels_, iteri)
    {
        const multiphaseInter::phaseModel& phasei = iteri()();

        auto iterk = iteri;

        for (++iterk; iterk != phaseModels_.cend(); ++iterk)
        {
            if (iteri()().name() != iterk()().name())
            {
                const multiphaseInter::phaseModel& phasek = iterk()();

                // Phase i and k
                const phasePairKey keyik
                (
                    phasei.name(),
                    phasek.name(),
                    false
                );

                if (interfacePorousModelTable_.found(keyik))
                {
                    autoPtr<porousModel>& interfacePtr =
                        interfacePorousModelTable_[keyik];

                    vector solidVelocity = get<vector>("solidVelocity"); //Edit to read from Subdict

                    Udiag += Vc*interfacePtr->S();
                    Usource += solidVelocity*(Vc*interfacePtr->S());
                }
            }
        }
    }
}

The same approach has been applied in the solver I was referring to (directChillFoam).


I hope this information proves to be useful for someone

[Alczem]

Hey!


I have also been looking into unconstrained solidification and melting (so a liquid could solidify and the resulting solid would keep flowing).


Considering your approach, does it mean that the solid phase remains the same shape but is still affected by the forces exerted by the liquid moving around it and does it still move accordingly? Can it sink because of gravity?


Thanks

[Lang]

Hello Alczem,

in this approach the solid phase maintains its shape and isn't influenced by external forces. Thus, it won't be affected by gravity and won't sink. The velocity field within the solid area is set to a predetermined value.
While it's theoretically possible to modify the code so that the forced velocity value varies over time, I don't believe this would address your specific problem. Perhaps someone else could provide more useful insights.

Best regards and good luck with your project,
Matthias