[palarcon]

Hello everybody

I'm working on topology optimization with OpenFOAM and I'm using and adapting the solver adjointShapeOptimizationFoam.
In order to do that, I need to compute and print the objective function used on that solver (which is not the same used on the paper which is based on).
The real objective function used in this solver is the total pressure drop, which mathematically is



which is the total pressure at the inlet minus the total pressure at the outlet.
I'm not sure about the implementation into my solver, mainly for the differential part .

My current implementation is

objective_inlet  = sum((p.boundaryField()[objPatchList[0]] + 0.5*magSqr(U.boundaryField()[objPatchList[0]])));
objective_outlet = sum((p.boundaryField()[objPatchList[1]] + 0.5*magSqr(U.boundaryField()[objPatchList[1]])));
objective = objective_outlet - objective_inlet;

where objPatchList[0] is the inlet and objPatchList[1] is the oulet.

My question is, should I multiply objective_inlet and objective_outlet by mesh.magSf()[objPatchList[X]] in order to account for the differential in the analytical expression?

PS: I know that I posted a similar question some time ago, but my knowledge of the problem was different at that time.

[FlyBob91]

You should multiply for the cell face area, sum all the contributes and finally divide for the area of the boundary. If you multiply only one time the area you change the dimension of your cost function.
If you need help I can post my code.

Best regards

[palarcon]

Hello FlyBob

Finally I did (a few days ago) something similar to what you suggest. Nevertheless, given that OpenFOAM coding per se is not my strong point, I would be really grateful if you could post your code.

Thank you in advance

Originally Posted by FlyBob91

You should multiply for the cell face area, sum all the contributes and finally divide for the area of the boundary. If you multiply only one time the area you change the dimension of your cost function.
If you need help I can post my code.

Best regards

[FlyBob91]

scalar jDissTotalOutlet(0);
scalar jDissTotalInlet(0);

label PatchIDInlet = mesh.boundaryMesh().findPatchID("inlet");
scalar areaInlet = gSum(mesh.magSf().boundaryField()[PatchIDInlet]);

label PatchIDOutlet = mesh.boundaryMesh().findPatchID("outlet");
scalar areaOutlet = gSum(mesh.magSf().boundaryField()[PatchIDOutlet]);


jDissTotalOutlet = 

sum(mesh.magSf().boundaryField()[PatchIDOutlet]* (p.boundaryField()[PatchIDOutlet] + 0.5 * rho.boundaryField()[PatchIDOutlet] * magSqr(U.boundaryField()[PatchIDOutlet]))) / areaOutlet;	


jDissTotalInlet = 
sum(mesh.magSf().boundaryField()[PatchIDInlet] * (p.boundaryField()[PatchIDInlet] + 0.5 * rho.boundaryField()[PatchIDInlet] * magSqr(U.boundaryField()[PatchIDInlet]))) / areaInlet;	

scalar diff = jDissTotalInlet - jDissTotalOutlet;

Info<<"Total Pressure at the inlet: "<<jDissTotalInlet<<endl;
Info<<"Total Pressure at the outlet: "<<jDissTotalOutlet<<endl;
Info<<"The total difference is: "<<diff<<endl;

Since my formulation is compressible, while I imagine your is incompressible, you don't need the density (rho) so you can delete that term.
Hope it helps

Best regards