Debug: Database connection successful
You are not logged in.
This run from zero to 4 was carried out in pursuit of the cause of incorrect thrust reported from initial tests of the NewNozzle case.
We have NO nozzle installed, but we were getting upwards of 30 force-tons of thrust ... OK ... we have an error in our configuration.
However, the ** nature ** of the error is not obvious, or it would have been spotted before the run was first launched.
The working theory in play is that a feature of OpenFOAM that provides automatic input of mass may be the source of the incorrect results.
The feature is nice... it allows the operator to specify a flow of mass into the model regardless of whatever else is going on. I've been using this feature from the beginning, since it was available and it definitely does what it is advertised to do. However, unbeknownst to me, this feature had the non-obvious side effect of adding energy to the situation without advertising. In the tail of the log of a run from 0 to 4, with a new report added, we can see that pressure at the intake is greater than the average pressure in the pipe. A possible explanation (not yet confirmed) is that the intake feature is quietly pushing new hydrogen into the heating pipe against the existing pressure, and in so doing helping to create thrust at the exit.
Another way of thinking about this is the energy budget... Our design calls for an energy budget of 40 MW. However, if the intake feature is doing what i think it is doing, it is adding some amount of energy to the budget in order to force hydrogen into a chamber that is already at a high level of pressure. Whatever that energy may be is not accounted for at present.
In a Real Universe engine, the energy required to pump propellant into a combustion chamber against the very high pressures already present comes from the propellant budget, except for the very few rockets that use electric (battery powered) pumps.
Here is the data from the end of a run from zero to 4 seconds, using the standard 2 kg/s feed and 40 MW of budgeted energy.
ExecutionTime = 11200.03 s ClockTime = 11207 s << 4 seconds in under 4 hours
volFieldValue p_avg write:
volAverage(heaterCells) of p = 720582.8 << Average pressure in all 9600 cellsvolFieldValue T_avg write:
volAverage(heaterCells) of T = 897.7824 << Average temperature in 9600 cellsvolFieldValue rho_mass write:
volIntegrate(heaterCells) of rho = 4.952744 << Total mass - most of initial load should be gone at 4 secondsvolFieldValue T_min_global write:
min() of T = 26.53625 at location (-239.5 0.175 2.292485e-17) in cell 8880 << Incoming gas is at 20 KelvinvolFieldValue T_max_global write:
max() of T = 2762.275 at location (4.5 0.195 0) in cell 9799 << Maximum is at exit after all heatingsurfaceFieldValue massFlow_in write:
sum("inlet") of phi = -1.999996 << Intake is 2 kg/s as plannedsurfaceFieldValue massFlow_out write:
sum("outlet") of phi = 1.773396 << Exit flow is below intake so mass is accumulating in pipesurfaceFieldValue p_inlet_avg write:
areaAverage("inlet") of p = 735063.4 << This is the new measurement.<< My interpretation: There is more pressure at the intake due to intake feature forcing gas in
End
ChatGPT5.2 has a recommendation for more information gathering.
inletForces
{
type forces;
libs ("libforces.so");// Only measure force on the inlet patch:
patches (inlet);// Field names (match your case)
pName p;
UName U;// Compressible: tell it which density field to use
rhoName rho;
rhoInf 1; // reference (won’t matter much here)// Optional: set a reference point (moments use it; forces don’t care)
CofR (0 0 0);// Output control:
writeControl writeTime;
writeInterval 1;
log true;
}
(th)
Offline
Like button can go here
A run from 4 seconds to 8 seconds just ended. This run included additional data gathering functions.
I don't know how to interpret the results but hope there is something useful here:
ExecutionTime = 11246.17 s ClockTime = 11269 s
volFieldValue p_avg write:
volAverage(heaterCells) of p = 665639.9volFieldValue T_avg write:
volAverage(heaterCells) of T = 757.2487volFieldValue rho_mass write:
volIntegrate(heaterCells) of rho = 4.465979volFieldValue T_min_global write:
min() of T = 26.83986 at location (-239.5 0.175 2.292485e-17) in cell 8880volFieldValue T_max_global write:
max() of T = 1789.563 at location (-0.5 -0.195 -2.162287e-17) in cell 239surfaceFieldValue massFlow_in write:
sum("inlet") of phi = -2surfaceFieldValue massFlow_out write:
sum("outlet") of phi = 1.893913surfaceFieldValue p_inlet_avg write:
areaAverage("inlet") of p = 676283forces inletForces write:
sum of forces:
pressure : (-27051.32 0 0)
viscous : (1.319819e-06 9.629071e-18 9.153922e-30)
porous : (0 0 0)
sum of moments:
pressure : (0 -4.582672e-15 1.051603e-11)
viscous : (3.828871e-24 2.240422e-25 -2.313687e-15)
porous : (0 0 0)End
(th)
Offline
Like button can go here
I was hoping that this had been read already but here it is About OpenFOAM
A solar thermal optical plane space vessel uses focused sunlight via large mirrors/concentrators, often with optical fibers, to heat a propellant (like hydrogen) to extreme temperatures, creating thrust, offering high efficiency (high specific impulse) for orbit changes, unlike solar sails that use light pressure, with concepts like Portal Space Systems' Supernova aiming for high maneuverability for defense/service, building on NASA's research for advanced solar propulsion systems.
Core Concepts
Solar Thermal Propulsion (STP): Captures solar energy with large mirrors (concentrators) and focuses it onto a heat exchanger/absorber, heating a propellant to generate thrust, achieving higher efficiency than chemical rockets.
Optical Waveguides: Some designs use optical fibers to transmit concentrated solar energy from the collector to the engine, allowing for more flexible system designs.
Optical Plane: Refers to the precise optical systems (mirrors, lenses) used to concentrate sunlight and potentially for imaging/sensing, requiring high-temperature materials and stability.
Key Components & Technologies
Concentrators: Large, often inflatable mirrors that focus sunlight.
Absorber/Heat Exchanger: Receives concentrated solar energy and heats the propellant.
Propellant: Typically hydrogen, heated to over 2500K, sometimes over 3000K with advanced materials.
Materials: Requires high-temperature refractory metals (tungsten, rhenium) or carbides to withstand extreme heat.
Optical Structures: Precision-built composite structures for stability in imaging and concentrating systems, used in telescopes and rovers.
Applications & Examples
Upper Stages: Ideal for moving payloads from Low Earth Orbit (LEO) to higher orbits or Earth escape.
Small Satellites: Enables orbit changes and long-duration missions with high delta-V.
Portal Supernova: A current project using solar thermal propulsion for high-delta-V, in-orbit refueling, and autonomous drone swarms.
NASA Research: Historically studied by NASA Marshall Space Flight Center (MSFC) and Sandia National Labs for advanced propulsion.
Contrast with Solar Sails
Solar Sails: Use direct photon pressure from sunlight on large sails, no propellant needed (e.g., IKAROS, LightSail-2).
Solar Thermal: Heats a propellant for rocket thrust; distinct from solar sails but both use solar energy for propulsion
Online
Like button can go here
Solar thermal propulsion (STP) utilizing a hydrogen propellant, enhanced by heat pipe technology and leveraging high-performance, regeneratively cooled engine techniques (similar to the Merlin engine's regenerative cooling), is a highly efficient propulsion concept for in-space transportation. This approach, often characterized as a Solar Powered Rocket with Impulsive Thermal Engine (SPRITE), can achieve specific impulse (\(I_{sp}\)) values of 830–1000+ seconds, roughly twice that of the best chemical engines. Key Components &
Technology Hydrogen Propellant: Used for its low molecular weight, providing superior \(I_{sp}\) (800-1000+ seconds) compared to other propellants.
Heat Pipe/Thermal Storage: High-temperature sodium wicking or phase-change materials (PCM) are used to store solar energy, enabling the system to "coast" and accumulate energy, then release it for high-thrust, intermittent burns.
Materials: Refractory metals such as tungsten or rhenium are required to handle high temperatures (>2500K) and prevent chemical degradation from hot hydrogen.
"Merlin" Capabilities (Regenerative Cooling): Similar to the SpaceX Merlin engine, which uses regenerative cooling to handle extreme temperatures, solar thermal engines use the cold hydrogen propellant to cool the engine chamber and nozzle before it is heated and expelled.
Performance Characteristics Thrust & Efficiency: While early designs produced around 8.9 N of thrust, advancements in high-temperature materials and concentrator technology are aimed at increasing performance.
Energy Collection: Uses deployable solar concentrators to gather solar energy over time, eliminating the need to have collectors pointed at the sun during all maneuvers.
Propellant Management: To handle cryogenic hydrogen, advanced zero-boil-off (ZBO) techniques are used to ensure long-term storability.
Operational Advantages High Specific Impulse (\(I_{sp}\)): Over 800s to 1000s, significantly higher than chemical (\(I_{sp}\) ~300-450s).
No Hazardous Materials: Unlike nuclear thermal propulsion (NTP), solar thermal systems are generally non-toxic and safer to operate.High-
Thrust Capability: The use of energy storage/impulsive thermal engines allows for high-thrust maneuvers suitable for orbital transfers (e.g., LEO to GEO in ~1 month).
Versatility: The system is ideal for in-space tugs, utilizing on-orbit refueling or potentially harvested water/hydrogen as propellant. Current research, such as the NASA-supported SPRITE project, is focused on reducing the total wet mass of these vehicles to below 200 kg while maximizing propellant efficiency
Solar thermal propulsion (STP) using hydrogen heat pipes represents a high-efficiency alternative to traditional chemical engines like SpaceX's Merlin. While STP systems excel in fuel efficiency (Isp), they operate at vastly different scales than the high-thrust Merlin engines used for heavy lifting.
1. Comparative Capabilities (STP vs. Merlin)
As of 2026, experimental STP systems are achieving performance metrics that surpass the efficiency of Merlin engines while trailing significantly in raw power.Specific Impulse (Efficiency): STP systems such as the SPRITE (Solar Powered Rocket with Impulsive Thermal Engine) achieve an Isp of ~830–860 seconds. Advanced designs using carbide materials aim for >1,000 seconds. In contrast, the Merlin 1D vacuum engine has an Isp of ~348 seconds.
Thrust-to-Weight Ratio: The Merlin 1D holds a world record thrust-to-weight ratio of >150. While STP systems like SPRITE are designed for high thrust-to-weight relative to other solar concepts, they generally produce low absolute thrust (e.g., 8.9 N) suitable for orbital transfers rather than Earth-to-orbit launches.
Mission Utility: STP is primarily targeted for upper-stage maneuvers, such as moving payloads from LEO to GEO or the Moon, where high efficiency over long durations is more valuable than the immediate raw thrust required by first-stage engines like the Merlin.2. Role of Hydrogen Heat Pipes in 2026
In 2026, hydrogen heat pipes are critical for thermal management in these high-temperature systems:Thermal Transport: Heat pipes function as "thermal superconductors," using phase-change fluids to move massive heat loads from solar concentrators to the engine's absorber cavity with minimal temperature drops.
System Enhancement: Recent developments show that integrating parabolic reflectors with advanced heat pipes can increase local temperatures by over 30%, significantly improving the melting performance of phase-change materials used for energy storage.
Propellant Management: Solar thermal rockets utilize hot hydrogen as a propellant because its low molecular weight allows for much higher exhaust velocities than the kerosene/oxygen mix used in Merlin engines.
3. Current Milestones (2025–2026)
Full-Scale Tests: In late 2025, companies like Portal Space Systems validated STP architectures at temperatures of 1,500°C (2,700°F) using electrical induction to simulate solar power, with full-scale demonstration missions scheduled for 2026.
Advanced Materials: Research in 2026 continues to focus on refractory metals (tungsten, rhenium) and high-temperature carbides to withstand the 2,500K+ environments required for efficient hydrogen heating
Online
Like button can go here
For all ... this topic is intended to serve as a tutorial for persons who are studying OpenFOAM for school, work or for personal projects.
We have attempted to report our experiences with OpenFOAM with as much detail as possible.
We have reported some successes, along with many failures.
A person should be able to read the topic from the top, and learn about various kinds of situations where OpenFOAM works well and others where it does not.
As of today, Wednesday 2026/01/21, we have a well debugged model that shows the original Merlin engine running with a 240 meter hydrogen heater in place of the original Merlin combustion chamber. We are working on a modification with a smaller nozzle throat, better suited for use with hydrogen as the propellan.
At the moment, we are debugging a case that is reporting ridiculous levels of thrust. The most likely cause for this misbehavior is a feature of OpenFOAM that generates artificial energy flow without the experimenter being aware of it. The feature is no doubt useful, but it is not helpful in every case.
For a practical engine, we will need to use some of the solar power (40 MW in the present case) to drive the pumps, so whatever is used to drive the pumps is not available to heat the gas.
(th)
Offline
Like button can go here
It is time to begin thinking about adding a nozzle.
Up to now, our model consists of cubes 10 centimeters on a side. We have 9600 cubes making up the heating pipe, and 200 making up the stub. ChatGPT5.2 reminded me that when we introduce curves we will need to greatly reduce the X interval. At the moment I am thinking of millimeter slices for curves, and centimeter slices for straight sections.
In this part of the model, I see two curves at a minimum. The first is the in-bending curve to make the transition from the heating tubes to the straight wall that leads to the throat. The throat is one continuous curve that begins at the slope from the heating tubes, and ends after an approximately 120 degree bend to join a straight section for the expansion wedge.
Guesses: 30 degrees bend to straight, 120 degrees at throat, 30 degrees straight for expansion wedge
Please note the terminology. We will not have an expansion bell, because we have a 40 centimeter wide heating tube feeding the nozzle. Rather than try to collect all the gas from the heater into a cylinder, we have chosen to keep the 40 centimeter width of the heater through the nozzle.
The gap at the throat will be something like 7 millimeters. The idea is to achieve the area for the throat equivalent to the area of a circle of 3 cm diameter.
GW Johnson has provided a drawing to show the design he would like to see. This location is available if we can obtain a copy of the image.
Please note that in the version I am planning for the NewNozzle case, I will attempt to provide a circular merge between the heat pipe wall and the 30 degree slope down to the throat. While GW indicates we don't have to worry about that, i do not want to see the turbulence in the OpenFOAM animation that square corners would cause. The consume machine cycles and thus extend run time.
(th)
Offline
Like button can go here
This post is about the latest run of NewNozzle from 8 seconds to 9 seconds.
ChatGPT5 is hot on the trail of our mysterious bug. It estimated the power needed to feed propellant into the heater at about 150 KW. I can understand that estimate, because the area of the pipe needed to feed hydrogen into the intake is small compared to that huge 40 centimeter by 10 centimeter back wall. All the thrust produced by the heater is pushing against that forward wall, so the propellant pumps have to push against that to inject gas into the chamber. However, since the area of the feed pipe is small compared to the total area of the forward wall, the pump will only be pushing against a fraction of the pressure on that forward wall. Since 150 KW is so much less than 40 MW I think we can ignore it for now.
here is one of the new functions we are adding to the case:
p_outlet_avg
{
type surfaceFieldValue;
libs ("libfieldFunctionObjects.so");
regionType patch;
name outlet;
operation areaAverage;
fields (p);
writeControl timeStep;
writeInterval 1;
log true;
}(th)
Offline
Like button can go here
We have results from the run from 8 to 9 seconds.
I don't know what the results mean but someone else may.
ExecutionTime = 2844.41 s ClockTime = 2847 s
volFieldValue p_avg write:
volAverage(heaterCells) of p = 669188.7volFieldValue T_avg write:
volAverage(heaterCells) of T = 738.8852volFieldValue rho_mass write:
volIntegrate(heaterCells) of rho = 4.551232volFieldValue T_min_global write:
min() of T = 26.8307 at location (-239.5 0.175 2.292485e-17) in cell 8880volFieldValue T_max_global write:
max() of T = 1741.948 at location (-0.5 -0.195 -2.162287e-17) in cell 239surfaceFieldValue massFlow_in write:
sum("inlet") of phi = -2surfaceFieldValue massFlow_out write:
sum("outlet") of phi = 1.945187surfaceFieldValue outlet_U_avg write:
areaAverage("outlet") of U = (443.816 -1.602165e-08 1.816319e-15)surfaceFieldValue p_outlet_avg write:
areaAverage("outlet") of p = 658078.8surfaceFieldValue p_inlet_avg write:
areaAverage("inlet") of p = 679982.6forces inletForces write:
sum of forces:
pressure : (-27199.3 0 0)
viscous : (1.313604e-06 3.83113e-18 9.104782e-30)
porous : (0 0 0)
sum of moments:
pressure : (0 -4.607741e-15 5.712764e-12)
viscous : (3.803829e-24 2.229902e-25 -9.209213e-16)
porous : (0 0 0)End
ChatGPT5.2 and i are on a quest to try to understand the invalid results we are seeing from a previously (apparently) reliable program "thrustCalc".
We moved a lot of the extendedMerlin case over to the NewNozzle case, and there may be further adjustments needed. Right now, thrustCalc is showing ridiculous numbers for thrust. I understand that physics equations predict how much thrust is possible given 2 kg/s of hydrogen and 40 MW of energy.
We have seen figures that match the predictions when we ran extendedMerlin. However, with NewNozzle we are seeing ridiculous numbers.
What is curious about this situation is that we do NOT have a nozzle installed yet. This was supposed to be a shake-down run with nothing but a heat pipe and an empty stub where the nozzle should be. I am looking forward to finding out what-the-heck happened.
Here is the plot generated from the data provided by the run of NewNozzle from 0 to 9 seconds. The plot does not reflect reality. On the other hand, the plot itself accurately reflects the numbers it was given.
(th)
Offline
Like button can go here
This is a detective story .... ChatGPT5.2 and I are following the clues to try to find the cause of the ridiculous results shown in the thrust plot.
We ** just ** learned that we had created the program some time ago. We both forgot about it. ChatGPT5 forgot because it is limited in the amount of memory it's managers allow it to retain. I forgot because I'm human.
Here is the likely cause of the incorrect results!
cat thrustCalc.C
// thrustCalc.C — OpenFOAM v12-compatible, standalone thrust post-processor
// Computes force on the FLUID across a patch:
// T = ∫_patch [ rho*(U·nHat)*U dA + (p - pAmb)*Sf ]
// Vehicle thrust = negative of this vector.
#include "argList.H"
#include "timeSelector.H"
#include "Time.H"
#include "fvMesh.H"
#include "volFields.H"
#include "surfaceFields.H"
#include "fvc.H"
#include "IOdictionary.H"
#include "polyBoundaryMesh.H"
#include "word.H"
#include "vector.H"
#include "OFstream.H"
#include "OSspecific.H" // mkDir
#include <fstream>
#include <iomanip>
using namespace Foam;
int main(int argc, char *argv[])
{
// Enable time-selection flags before args are created by setRootCase.H
timeSelector::addOptions();
// Nice help text with --help
argList::addNote
(
"Compute thrust on a patch:\n"
" T = ∫_patch [ rho*(U·nHat)*U dA + (p - pAmb)*Sf ]\n"
"Force reported is on the FLUID; vehicle thrust is its negative."
);
#include "setRootCase.H"
#include "createTime.H"
instantList timeDirs = timeSelector::select0(runTime, args);
#include "createMesh.H"
// Read settings from system/thrustCalcDict
IOdictionary dict
(
IOobject
(
"thrustCalcDict",
runTime.system(),
runTime,
IOobject::MUST_READ,
IOobject::NO_WRITE
)
);
word patchName("outlet");
dict.readIfPresent("patch", patchName);
word inletName(""); // leave absent to ignore
dict.readIfPresent("inletPatch", inletName);
scalar pAmb(0.0);
dict.readIfPresent("pAmb", pAmb);
// Output path
fileName outDir = runTime.path()/"postProcessing"/"thrustCalc"/patchName;
mkDir(outDir);
const fileName outFile = outDir/"thrust.dat";
// Helper: find patch index by name (v12-safe)
auto findPatch = [&](const word& name) -> label
{
if (!name.size()) return -1;
const fvBoundaryMesh& bnd = mesh.boundary();
forAll(bnd, i)
{
if (bnd[i].name() == name) return i;
}
return -1;
};
forAll(timeDirs, ti)
{
runTime.setTime(timeDirs[ti], ti);
mesh.readUpdate();
// Exact time directory name (no rounding)
const word curTName = timeDirs[ti].name();
// Read fields at this time
volVectorField U
(
IOobject("U", curTName, mesh, IOobject::MUST_READ, IOobject::NO_WRITE),
mesh
);
volScalarField p
(
IOobject("p", curTName, mesh, IOobject::MUST_READ, IOobject::NO_WRITE),
mesh
);
volScalarField rho
(
IOobject("rho", curTName, mesh, IOobject::MUST_READ, IOobject::NO_WRITE),
mesh
);
// Compute contribution on a patch (fields are in scope here)
auto computePatch = [&](const word& name, bool subtract) -> vector
{
if (!name.size()) return vector::zero; // silent ignore
vector result(vector::zero);
const label pid = findPatch(name);
if (pid < 0)
{
Info<< "thrustCalc: patch '"<< name <<"' not found (skipping)." << nl;
return result;
}
const vectorField& Sf = mesh.Sf().boundaryField()[pid];
const fvPatchVectorField& Ub = U.boundaryField()[pid];
const fvPatchScalarField& pb = p.boundaryField()[pid];
const fvPatchScalarField& rb = rho.boundaryField()[pid];
const label N = Sf.size();
for (label i=0; i<N; ++i)
{
const scalar magSf = mag(Sf[i]);
const vector nHat = (magSf > VSMALL ? Sf[i]/magSf : vector::zero);
const scalar Un = (Ub[i] & nHat);
const vector mom = rb[i]*Un*Ub[i]*magSf; // rho * Un * U * dA
const vector pTerm = (pb[i] - pAmb)*Sf[i]; // (p - pAmb) * Sf
result += (mom + pTerm);
}
return subtract ? -result : result;
};
const vector thrust = computePatch(patchName, /*subtract*/false)
+ computePatch(inletName, /*subtract*/true);
const scalar t = runTime.value();
const scalar Fmag = mag(thrust);
// Decide if we need a header (file does not exist or empty)
bool needHeader = true;
{
std::ifstream chk(outFile.c_str());
needHeader = !chk.good() || (chk.peek() == std::ifstream::traits_type::eof());
}
// Append one row (with header if needed)
std::ofstream os(outFile.c_str(), std::ios_base::app);
if (needHeader)
{
os << "# timeDir t[s] Fx[N] Fy[N] Fz[N] |F|[N]\n";
}
os.setf(std::ios_base::scientific);
os << curTName << " "
<< std::setprecision(10) << t << " "
<< thrust.x() << " " << thrust.y() << " " << thrust.z() << " "
<< Fmag << "\n";
os.close();
Info<< "thrustCalc: t="<< t
<< " Thrust_on_fluid [N] = " << thrust
<< " |F|=" << Fmag
<< " (vehicle reaction = -value above)" << nl;
}
Info<< "thrustCalc: wrote " << outFile << nl;
return 0;
}Update Friday 2026/01/23 ... I gave the program to ChatGPT5.2 to examine, and it came back with a suggestion that the error we are seeing might involve a sign going the wrong way in one calculation. However, the more interesting observation is that the program is reporting pressure on the end caps, and NOT thrust as it was supposed to be doing. I am shaky on this so am not sure, but I ** think ** thrust should be discoverable if we know the amount of mass we are sending out the exit port, and the velocity of that mass. I'm pressing ChatGPT5.2 to look into that question, along with numerous others.
Note: I may be a customer with the longest continuous project on record with ChatGPT... I don't have any way of ** knowing ** that. What I **know** for an absolute fact is that my queries are now taking literally hours to process, and when the do come back, they are marked with times such as "57 seconds" I think that 57 seconds may be machine time. In other words, it is possible that I am really getting my money's worth, but the trade is time between results.
(th)
Offline
Like button can go here
This is a follow up to #284
It took a long time (clock time) to work through the problem reported in #284, but ChatGPT5.2 and I decided to re-compile thrustCalc.C with just one change. We changed the line that works with the pressure readings at the end caps. The first program assumed that the pressure reading at the intake is a positive value, but in our case it is negative. The result was that the thrust plot added the pressures instead of subtracting them. That led to the ridiculous plot that I've been chasing ever since. I just recompiled the program and it is ready to run just as soon as I send this report.
The thrust from our heat pipe with no nozzle should be very low... ChatGPT5.2 thinks we should see something like 800 Newtons.
We spent a ** lot ** of time discussing choked flows. GW Johnson understands this and despite ChatGPT5 and Gemini's best efforts, I ** still ** don't understand it, but I accept that whatever this is, it is real and measurable and every ** real ** rocket scientist believes it.
I don't understand it because, unlike an engineer, I want to understand the physics behind the observations. An engineer doesn't worry about the physics. A ** real ** engineer just takes the formula and runs with it.
Here is the corrected plot of thrust for the NewNozzle case:
We have a 240 meter long ventilation duct that is 40 centimeters wide and 10 centimeters high.
It is filled with hydrogen at 1 bar and 21 Kelvin. When we begin the run, we open the exit to near vacuum and we begin heating with 40 MW while we ramp up a flow of hydrogen at the intake from 0 kg/s to 2 kg/s over 2 seconds. What I think the plot now shows is the hydrogen rushing out the exit while pushing gently on the hydrogen in the pipe so that a very small residual force is expressed at the intake.
I am looking forward to hearing feedback on this new set of results.
I should emphasize that this was NOT caused by a bug in the program thrustCalc.C
This was a switch setting under the control of the operator. The fact I was unaware of this switch and the incorrect setting does not let me off the hook. What I ** can ** claim is sufficient understanding to know that earlier results were wrong, and (more importantly) the persistence to keep digging until ChatGPT5.2 and I finally arrived at a consensus on the nature of the problem (a) and (b) the solution.
(th)
Offline
Like button can go here
This post contains a closing comment from the project that just ended. A "project" is a series of exchanges between a customer and ChatGPT, as distinct from a Chat, which normally runs for a short time such as a day or two at most. A project can extend for months. In this case, we have finished development of extendedMerlin and successfully transitioned to NewNozzle (with some difficulty). It is time for a new project devoted to NewNozzle.
I am hoping that starting a new project will allow response time to return to normal. Right now, it takes hours of clock time between query/responses.
Here is the concluding comment by ChatGPT5.2 ...
tahanson43206 — this is a huge milestone. And yes: what you’re seeing now is exactly what “winning through the bramble patch” looks like.
1) Is the problem in your bash/awk/gnuplot pipeline?
No — your downstream pipeline is fine.
thrustCalc writes thrust.dat
your awk just flips signs for the components (vehicle reaction) and computes convenient units
gnuplot plots column 7
So if the curve was “ridiculous” before, it was because thrust.dat contained a quantity that was being misinterpreted as thrust (it was effectively “end-cap pressure force / control-surface force”, not the net propulsive thrust you intended). Now that you fixed the arithmetic, the CSV is reflecting the corrected quantity.
2) Does “one-line arithmetic fix” make sense?
Yes — and your new thrust_series.csv is the proof. The original symptom (tens of kN in a straight open duct with no nozzle) is exactly what you get when the code accidentally adds pressure contributions that should cancel (or should be subtracted with the correct sign convention).
3) Why is the outlet showing pressure at all, if it’s “open to vacuum”?
Because your outlet condition is not “hard vacuum pressure = 0 Pa.” It’s a numerically stable, non-reflecting outlet that lets the interior solution determine the boundary pressure most of the time.
In your case the duct is running subsonic and the entire pipe is sitting at a high static pressure (hundreds of kPa). For subsonic outflow, the interior can “communicate upstream,” and the outlet patch will happily sit at almost the same pressure as the duct.
Also: the value uniform 20000; you set long ago is a fallback used only in certain backflow / initialization situations. It is not a clamp. So it’s completely normal to measure ~600–700 kPa on the outlet patch even though the external world is “vacuum” in your conceptual vehicle.
4) “Why are we using pressure at all to compute thrust?”
Because the control-volume thrust balance has two parts:
momentum flux (mass flow and velocity)
pressure forces acting on the control surface
A common “rocket thrust” identity is:
F = mDot * Vexit + (Pexit − Pamb) * Aexit
That second term is the pressure term. It’s not optional.
If the jet is perfectly expanded to ambient, then Pexit ≈ Pamb, so that term goes away. But in your straight duct case, your “exit” is not a properly expanded nozzle flow — so the pressure term can be large on paper.What you discovered (the hard way) is the subtle but crucial point:
The pressure term must be included correctly so the inlet and outlet pressure contributions cancel appropriately in a duct-like control volume.
If you accidentally “double count” it, you can fabricate enormous forces that aren’t real propulsion.
5) So what is thrustCalc reporting now — “true thrust” or something else?
With your corrected arithmetic, it is now much closer to what it should be: a net streamwise force balance across the chosen control surfaces (outlet, optionally corrected by inlet).
That’s why you now see small values (tens of N) at late time — that’s consistent with “a duct on frictionless ice” where the net propulsive force is basically the net momentum change of the flow plus any non-canceling pressure term. In a straight constant-area duct with similar p at inlet/outlet, the pressure terms nearly cancel, leaving a small net.
6) Can we “fix thrustCalc.C to give us true momentum thrust”?
Yes — and it’s actually easy conceptually:
If you want momentum-only thrust, compute only:
Fm = ∫ rho(U·nHat)U dA over the outlet, and (optionally) subtract the same at the inlet.
That will give you the “mDot*V” part directly.
But for your team’s purposes, I’d recommend outputting both every time step:
F_momentum (outlet − inlet)
F_pressure (outlet − inlet), using (p − pAmb)*Sf
F_total = F_momentum + F_pressure
Then nobody has to guess what’s being plotted.
7) Should we “start over” and write a new program?
We can, but you don’t need to. Your current program is already 95% of the way there — it just needs to write out the components separately (momentum vs pressure), and to be explicit in the header/labels.
One practical next step (one-variable-at-a-time, as you like)
Add two extra lines of reporting from the solver (not thrustCalc), so we can cross-check thrustCalc numerically:
sum(outlet) of phi (you already have)
areaAverage(outlet) of U (you already have)
and add areaAverage(outlet) of rhoThen you can do a back-of-envelope momentum thrust right in the shell:
Fm ≈ mDot_out * Ux_out
and compare that to thrustCalc’s momentum-only output (once you split it).
If you want me to sanity-check a specific timestamp numerically (using your measured mDot_out, Ux_out, p_out, and exit area), paste just one line of those four values and I’ll compute the expected momentum term and pressure term with your geometry.
The project we are closing today was called "OpenFOAM run Analysis"
Here is the corrected plot of thrust from the "NewNozzle" with no nozzle.
(th)
Offline
Like button can go here
This topic started in April of 2025. It was intended to provide an educational opportunity for anyone who might be interested in OpenFOAM, and secondarily, in a project to bring an inherited mesh of the SpaceX Merlin engine back to life. The topic includes considerable detail.
We close the phase of the project that led to creation of a brand new mesh for the Optical Plane vessel. The project will remain on-the-shelf until such time as there is an opportunity to resume.
In the mean time, other projects have been on hold. If there is a reader who is not a member of the forum, who would like to ask a question or make a contribution (such as a new nozzle mesh) please contact the forum using the address given in the Recruiting topic. To see the Recruiting topic, click on the Active button on the main menu.
(th)
Offline
Like button can go here