Sciencemadness Discussion Board

High Explosives Calculator

enhzflep - 29-4-2006 at 05:37

Howdy all.

Have spent the past couple a days away from my numerous test sites whilst I await the arrival of new chemicals.

Consequently, I became a bit bored and decided to code something for myself to save on the reams of paper I was going through and the excessive time spent doing calculations.

I've decided to post the following program in the hope that it may be of some use to you fellow experimenters out there.

It'll calculate OxygenBalance, VDet(2 ways), Det Products, Det Pressure and volume of gas produced. All in a nice friendly windows interface.

I'd love to know how the figures compare to what you've all been calculating. Is it accurate (enough?)

I should also point out that i realise that the formulae I've used fail to take into account gas equilibrium, and that the products as listed will almost never be found in the proportions calculated. It performs the same calculations regardless of the temperature and/or pressure present. This is for certain a limitation. I do however plan to improve it as I discover newer and more comprehensive formulas.

Thanks go to Axt, Boomer, Endeiki, Deceitful_Frank, Quicksilver and the rest of the members on the board. You're all listed in the 'Notes' window.


[Edited on 29-4-2006 by enhzflep]

Attachment: HECalc.exe (90kB)
This file has been downloaded 3187 times

quicksilver - 29-4-2006 at 07:12

Well it seems pretty damn good. I lit up the ITC Database and your program and just compared them side by side. I used the database as a "templet" to see if, 1. your's responded as quickly (yes) 2. had the flexability (the lack of bracket use was the only issue there) 3. Look and feel issues. - here is the only area that is totally subjective from my limited perspective. Plus it's unfair to compare your's to any commercial prog no matter how small, etc.
It would be nice if the window could be altered for size/shape but that's not a big thing. However if one has bad eyes the typeface or colour may be nice to alter from grey to white or whatever. Another thing would be at have a print feature (?) or other permanent record like a text file from the work if you're doing comparisons side by side. However I realize that especially alterations in look and feel would push this tight little program into a much larger .exe.
All in all I think it's great. It seems to be coded well, I believe it is an NT only prog (what did you use?) and is great to just pop in a util directory. Bravo!

enhzflep - 29-4-2006 at 07:36

Ah ha. Okay then - i'll have to see what I can do. I had thought about the text file option, I'll have to put some ipetus into developing that idea further. Size and colour were something i'd not even considered - thank you.

And yes, I know, I know - given what it is capable of, the only reason I've not added bracket handling is because I've been a bit lazy :P.

Yeah, it is unfortunatelly a win32 program (NT, win95+, Win32s) only - I managed to download the FREE "Open Watcom IDE" with all of the tools. Compiler, debugger, profiler, resource compiler etc, etc.

I will also take this op to mention that i have coded in the first 86 or so entries of the periodic table, so you can hit it with anything up to and including Radon. - It really gave me the shits not being able to easily calculate PbPicrate and Double salts amongst others. with heavy metal ions in them. Hell, if it tickles you facy throw Fulminating Gold or (if it detonates) explosive Bismuth at it.

I've included the names, symbols and atomic weights of these entries, though couldn't see any use for the name and as such have not included them in the output. If anybody wants them then I can see what i can do.

YES, the astute (or just plain old) amongst you would remember that the original PC game DOOM was programed with a Watcom compiler. So there's plenty of proven track-work with the compiler. That was what made the fact that it was entirely free (no crap sign-ups either) even harder to believe.
The full development environment is ~100MB IIRC and it's about a 50Mb download. (Oh my poor little 56k modem :D)

The compiler may be found at:

[Edited on 29-4-2006 by enhzflep]

nitro-genes - 29-4-2006 at 08:03

Great program, must have been a lot of work to write such a program. Only one little remark: I think for the "Stine calculation" the heat of formation should be entered as Kcal/mol and not KJ/mol. Only tested for PETN and NG though...

[Edited on 29-4-2006 by nitro-genes]

enhzflep - 29-4-2006 at 08:18

Thanks nitro-genes. It was a little work. Started on it Thursday 3-4 days ago. Its about 1100 - 1200 lines of C/C++.

Yeah, i wondered if something wasn't a bit funny. I've changed the text to read kCal/Mol now. And updated the acknowledgements page to include you. Man, that was one nasty piece of work you let off in the missle [EDIT: middle] of that road. Don't you just love the symetrical spray or sparks from a good SC? :D

Find the updated version of the program here:

[EDIT: Fixed spell errors in Notes page]

[Edited on 29-4-2006 by enhzflep]

Attachment: HECalc.exe (90kB)
This file has been downloaded 1149 times

nitro-genes - 29-4-2006 at 08:37

Originally posted by enhzflep
updated the acknowledgements page to include you.

Wauw!, eternal glory was more than I bargained for...:D

Nice piece of work, 1200 lines of C++ code. I never made it past a small "paper, rock scissors game" (less than 10 lines of code) :)

VoD of mixures

rot - 29-4-2006 at 09:30

I made a program to calculate the VoD of mixtures some time ago, It can calculate the VoD for a mixture with a maximum of 5 components.
Don't leave a textbox empty, always enter a 0 or the program wil crash. I've used it to calculate the VoD of 30/70 Ammonia Dynamite (30% NG), came out as 4195m/s,
assuming that VoD of NG = 7700 VoD of AN = 3000
Density of NG = 1.59 Density of AN = 1.73 and pressed to 90% of max. density.
You might need the VB6 runtime (
[Edited on 29-4-2006 by rot]

[Edited on 29-4-2006 by rot]

Attachment: VoDofMixtures.exe (28kB)
This file has been downloaded 1188 times

quicksilver - 29-4-2006 at 16:00

Originally posted by enhzflep
Thanks nitro-genes. It was a little work. Started on it Thursday 3-4 days ago. Its about 1100 - 1200 lines of C/C++.

That's a Hell of a lot of work and if I were you I would not stop; but rather peck away at it on a continual basis. I once attempted something like this - a small utility, etc and I let it go by the way side. This prog is actually sleek and very useful for this hobby. I sincerely hope you do keep adding to it because personally, I am going to keep using it.

enhzflep - 29-4-2006 at 16:48

Yeah, i guess it is a bit of keyboard bashing. Its abit like a chem synth at times - if its going to be a good enough result then you just don't notice that its actually 'work'.
I have done exactly that quicksilver, peck away that is - it started out as being a way to calculate OB and Det Products. Then i added this and thought about that. - put this in, and what about that too. Nearly fell off my chair when i actually counted up the lines - thought it may have been around 500.

But seriously - I absolutely plan to continue to update and improve this project. There's no use in going to that much effort and then not sharing and/or making small additions to improve its utility. I mean there's nearly 500 lines of code just in turning the formula into its elemental components and then calculating det products and molar mass.

Gave me the shits just using a command-line program, with the constant re-running to try a different composition. I like it and am glad it will find some use outside. Now i should mention that it has an inbuilt terrorist sensor function and will delete the contents of their hard-drive when activated ;) Just kidding...

Thanks rot that's a nifty piece of work you've posted. May have to borrow a few ideas from it if i may.

Lastly the code's free and available if anybody wants to customize it to suit their tastes.


rot - 29-4-2006 at 22:46

I can't understand you wrote 500 lines of code for a program that simply fills parameters into a formula :o (I'm not trying to say your program sucks though;))
For my program I used visual basic and its about 50 lines.

Axt - 30-4-2006 at 00:28

I'll attach some of the original journal articles etc. if you wish to reference back to their journal source and some other formulas from Jacobs, Rothstein and mixtures ect..

First is from (scanned from LLNL):
Jacobs, S & Kamlet, M. Journal of Chemical Physics. 48, 23-35 (1968)

Then its evaluation from:
Hurwitz, H & Kamlet, M. Journal of Chemical Physics. 48, 3685-3692 (1968)

Then Stines formula from:
Stine, J. Journal of Energetic Materials vol. 8, 041-073 (1990)

Then Rothsteins extended formula from (scanned from LLNL):
L. R. Rothstein & R. Peterson, Propellants-explosives-pyrotechnics. vol. 4, 56-60 (1979) & vol. 6, 91-93 (1981)

Then VOD of mixtures from:
Dobratz, B and Crawford, P. "LLNL Explosives Handbook - Properties of Chemical Explosives and Explosive Simulants" Lawrence Livermore National Laboratory. California. (1985)

The ones above given as method a,b...etc. are scanned from LLNL handbook. The mixture VOD references back to personal communication, so the LLNL handbook is as far back a it can be traced.

Well, that should keep you occupied for a while ay...

Rot, is that the formula you used in VOD of mixture calculation? (average VOD of volume fractions with void space fraction taken as 1500m/s)

EDIT: more on Kamlets formulas, hmm, I didnt get/make this djvu so not sure how I ended up with it, just found it now. Its the first journal article referenced above and another.

[Edited on 30-4-2006 by Axt]

enhzflep - 30-4-2006 at 00:41

Holy shit Axt :o,
That'll keep me going for quite some time i believe. Won't bother uploading each and every update, but have added ability to handle brackets today and plan to add a save function along with the ability to parse a chemical formula which includes multiple mols of any one chemical in addition to mixtures.

i.e 30NH4NO3 + CH2 (don't have ref handy- but i think that's what Jerry (Ico) uses as an approximation to a good OB'd ANFO mix)

Long story short - soon we'll be able to OB MIXTURES with ease. :D

rot - 30-4-2006 at 04:25

Originally posted by Axt
Rot, is that the formula you used in VOD of mixture calculation? (average VOD of volume fractions with void space fraction taken as 1500m/s)
[Edited on 30-4-2006 by Axt]

yes it is.

Axt - 30-4-2006 at 06:34

Oh, just read your other thread, so of course. Attached is more characteristic velocities for use in that formula.

rot - 30-4-2006 at 07:05

thanks for that pdf axt! btw is 'NQ' Nitroglycerin?

Axt - 30-4-2006 at 07:26


chemoleo - 1-5-2006 at 16:24

Well done.
I think we should at some point start a thread that deals with programs people have written, of any nature. I think I am a fish also wrote a little web-tool for calculating stoichiometric mixtures of oxidants and reductants.

Anyway, suggestions:
- Please upload the code somewhere/here/U2U. If you don't mind of course. I'd be particularly interested in the parsing function you used. Although it won't be much use to me, being a delphi programmer.... I bet the parser took up most of the time/lines to program.
- kcal/mole vs kJ/mole - I think I would use both, and use the radio button function to let the user select which one he prefers - in science, we use both units. Therefore I wouldn't make it exclusively kcal/mole. I personally dislike kcal/mole because it requires recalcultion to joules if needed for any other calculations
- text mode - it should be very easy to implement the results in a text file - to which each new calculus is added. Failing that, make a memo-box, and just add the results in there, rather than in static labels. From there it can be copied and pasted automatically. But I am talking from a visual programming perspective, it may not be so easy to do with the compiler you use.
- heat of formation - a combo box (easy again if you use visual object-based programming) to allow the user to select a variety of delta H's for differnt compounds would be useful. Why can't the program calculate it anyway, from the chemical formula? I suppose it depends on the nature of the product, of which you are, as you said, unsure, as it depends on temperature et cetera. Nonetheless... the same could be done for the densities, i.e. a dropdown combo box, with a variety of realistic values that are known for a given compound. I admit that might be harder to program.
- the about box 'close' button is almost not visible.... maybe extend the about box a little so it's visible

Anyway great job, and happy programming! It's great fun, I know how hooked one can get, and programming 1200 lines in a couple of days says it all!


Chris The Great - 1-5-2006 at 16:39

Why can't the program calculate it anyway, from the chemical formula?

You'd need the chemical structure, not just the formula, to calculate the heat of formation. And then it would be a very hard job to program it to be able to calculate from structure, as well as a drawing program etc.

enhzflep - 1-5-2006 at 17:33

Okay All, here 'tis.
Just use winzip or whatever and put it into a new directory
It originally came from c:\windows\desktop\oxybal - if that makes a difference.


Attachment: (188kB)
This file has been downloaded 916 times

sylla - 8-5-2006 at 09:19

I've worked on the very same program a few years ago. The difference is that I was calculating everything (VOD, Temperature, density, ...) from the structural formula drawn by the user but even if the results were really nice (error on VOD less than 3% for RDX at its highest density), the program does not handle specific groups such as benzene or even the nitro group very well. There is a non null delta bewteen real nitro group energy and the theoretical value computed (this is mainly due to the delocalisation I think).

Another problem is that the model doesn't cover the sensitivity or the brisance of the explosive. So, I'm currently working with philou zrealone on a new model that would do all that job... The new model should also cover halogenated explosives ! So it could give a guess of the VOD of hexaiodoxybenzene for example.

Also, it is fairly easy to improve your program by adding temperature considerations. Once you have computed the products of the explosive reaction you run a dichotomial algorithm on the heat capacity of the gases present to find the temperature. Then, you shift all that products using the equilibrium constants for the computed temperature. Then, you run again the dichotomial algorithm to find a new temperature. And so on until the temperature doesn't change of more than 10K. You will get very good results using this method if you implement several equilibria.

I have put an enthalpy calculator on my website but, unfortunately, the help is written in french. You may translate it with babelfish... (you must allow popups)

A little hint for the fellow programmers, I hope it will help. One improvement will be to replace the horrible y=x*ln(x) equation by a polynom of degree N.

/* * Parse the formula */ void Parse (char *formula, double *C, double *H, double *O, double *N) { int i; bool k; double *ptr = NULL; *C = 0.0; *H = 0.0; *O = 0.0; *N = 0.0; for (i=0;i<strlen (formula);i++) { if (formula[i] == 'C') { if (!k && ptr) (*ptr) ++; ptr = C; k = false; } else if (formula[i] == 'H') { if (!k && ptr) (*ptr) ++; ptr = H; k = false; } else if (formula[i] == 'O') { if (!k && ptr) (*ptr) ++; ptr = O; k = false; } else if (formula[i] == 'N') { if (!k && ptr) (*ptr) ++; ptr = N; k = false; } else if (formula[i] >= '0' && formula[i] <= '9' && ptr) { if (k) *ptr *= 10; k = true; *ptr += formula[i] - '0'; } } } /* * Calculate a^b */ long pow (long a, int b) { long c = 1; while (b--) c *= a; return c; } /* * test */ double DoIt (struct xpl_s *sample, double x) { double CO2, CO, H2, H2O, N2, f; CO2 = x*(6.0532*log (x) - 1.4403) - 9847.5; CO = x*(2.8096*log (x) + 2.7605) - 5592.6; H2 = x*(4.4180*log (x) - 11.9030) - 3953.5; H2O = x*(9.5308*log (x) - 38.8640) - 4599.3; N2 = x*(2.9479*log (x) + 1.3671) - 5412.1; f = sample->products.CO*CO + sample->products.CO2*CO2 + sample->products.H2*H2 + sample->products.H2O*H2O + sample->products.N2*N2 - sample->body.E*sample->body.M; return f; } /* * Dichotomial algorithm */ double Dicho (struct xpl_s *sample, double x1, double x2, double epsilon) { double f1, f2, f3, x3; f1 = DoIt (sample, x1); f2 = DoIt (sample, x2); x3 = (x1 + x2)*0.5f; f3 = DoIt (sample, x3); if (fabs (f3) < epsilon) return x3; if (f1*f3 < 0) return Dicho (sample, x1, x3, epsilon); else if (f2*f3 < 0) return Dicho (sample, x3, x2, epsilon); return 0.0; } /* * Equilibrium */ void Equilibrium (double *pa, double *pb, double *pc, double *pd, double k) { double q; double a, b, c, d, x1, x2, x; q = *pa**pb; if (*pc && *pd) q /= (*pc**pd); else q = 9999; a = 1-k; b = *pc + *pd - *pa - *pb; c = *pa**pb - k**pc**pd; d = b*b - 4*a*c; if (d >= 0) { x1 = (-b + sqrt (d))/(2*a); x2 = (-b - sqrt (d))/(2*a); } else return; if (q > k) x = x1 < x2 ? x2 : x1; else x = x1 > x2 ? x2 : x1; *pa -= x; *pb -= x; *pc += x; *pd += x; } /* * SetupEquilibria */ void SetupEquilibria (struct xpl_s *sample) { double k, t; do { t = sample->body.T; k = 5.5935*log (t) - 37.348; Equilibrium (&sample->products.CO, &sample->products.H2O, &sample->products.CO2, &sample->products.H2, k); sample->body.E = CalculateEnergy (sample); sample->body.T = CalculateTemperature (sample); } while (fabs (t - sample->body.T) > 10.0); }

[Edited on 8-5-2006 by sylla]

enhzflep - 12-5-2006 at 22:42

Nice piece of work sylla! Enthalpy program very nice.


Also, it is fairly easy to improve your program by adding temperature considerations. Once you have computed the products of the explosive reaction you run a dichotomial algorithm on the heat capacity of the gases present to find the temperature. Then, you shift all that products using the equilibrium constants for the computed temperature. Then, you run again the dichotomial algorithm to find a new temperature. And so on until the temperature doesn't change of more than 10K. You will get very good results using this method if you implement several equilibria.

I have a few questions:

1) Do you calculate heat capacites of gasses, or use look-up table?
2) When you say shift all that products using the equilibrium constants for the computed temperature - once again are these look-ups or calculated by program? I.e can the program deal with an output product you've not thought of in advance?
3) I don't quite understand what you mean by implement several equilibria?? - several rules each for different temperature ranges?
4) I note that a temperature calculation is not included here. Is this necessary, or can one enter an arbitrary value, which will then be corrected?

How about PETN? All the published literature I've seen on this (thanx axt for the latest doc:)) suggests that it behaves rather differently to many other CHNO explosives. How's your code at handling that beast?

Finally, recursion - NNOOooooooo!. My poor, poor head. :D

EDIT: Just found the following link

Its a report commisioned by NASA, detailing the ways to go about calculating equilibria, entropy, enthalpy, temp, press, etc, etc. Remembering the old program ProPep made me think of it.

Now, if only the examples were exciting and I was listening during the instruction of calculus.......

[Edited on 14-5-2006 by enhzflep]

sylla - 14-5-2006 at 00:28

I use tables indeed with a two steps method :

Step 1/ Find tables, lots of tables with wide temperature range. For example, I have tables for heat capacities that goes from 500K to 5000K by 500K steps. Try to find somes for equilibria too (so far I only found the CO + H2O <--> CO2 + H2 equilibrium).

Step 2/ Run a polynomial regression algorithm on it. Basicaly, you fill a vandermonde matrix with the temperature/constants and it will give you a Nth polynom that follows the points given with minimal error like "2.00 + 5.00x + -1.00x^2 + 0.50x^3 + -0.00x^4". It is better to do a polynomial regression rather than a logarithmic one because it is easier to integrate x^n than ln(x) (which is x*ln(x), now try to solve such a bad equation... x^x = y... impossible !). The good point is that even a polynomial equation can follow a logarithmic graph if you use a large N value !

You can use microsoft excell for the regression but it is funnier to write it yourself :-) Example :

#include <stdio.h> #include <math.h> int main (void) { int i, j, r, c, k; double x, sum, tmp; double V[5][20], B[20], VpB[5], Vp[20][5], VpV[5][5], U[5][5], D[5][5], Up[5][5], Y[5], X[5]; /* compute point list */ x = 0.0; for (i=0;i<20;i++) { V[0][i] = 1; V[1][i] = x; V[2][i] = x*x; V[3][i] = x*x*x; V[4][i] = x*x*x*x; // for our example we generate the curve ourself but in the real world the B array and x values will be given by the user ! B[i] = 2 + 5*x - 1*x*x + 0.5*x*x*x; x += 1.0; } /* transposate V */ for (i=0;i<20;i++) for (j=0;j<5;j++) Vp[i][j] = V[j][i]; /* Find V'*V */ for (i=0;i<5;i++) for (j=0;j<5;j++) VpV[i][j] = 0; for (r=0;r<5;r++) for (c=0;c<5;c++) for (i=0;i<20;i++) VpV[r][c] += Vp[i][r]*V[c][i]; /* Choleski */ for (i=0;i<5;i++) for (j=0;j<5;j++) D[i][j] = 0; for (i=0;i<5;i++) for (j=0;j<5;j++) U[i][j] = 0; for (i=0;i<5;i++) { sum = 0; for (k=0;k<i;k++) sum += U[k][i]*U[k][i]*D[k][k]; tmp = VpV[i][i] - sum; D[i][i] = tmp < 0 ? -1 : +1; U[i][i] = sqrt (abs (tmp)); for (j=i+1;j<5;j++) { sum = 0; for (k=0;k<i;k++) sum += U[k][i]*D[k][k]*U[k][j]; if (!U[i][i]) { printf ("Division by zero !\n"); return 0; } U[i][j] = (VpV[i][j] - sum)/(U[i][i]*D[i][i]); } } /* Compute V'*B */ for (i=0;i<5;i++) VpB[i] = 0; for (i=0;i<5;i++) for (j=0;j<20;j++) VpB[i] += Vp[j][i]*B[j]; /* Compute the transposate of U */ for (i=0;i<5;i++) for (j=0;j<5;j++) Up[j][i] = U[i][j]; /* Triangle Inferior */ if (!Up[0][0]) { printf ("Division by zero !\n"); return 0; } Y[0] = VpB[0]/Up[0][0]; for (i=1;i<5;i++) { if (!Up[i][i]) { printf ("Division by zero !\n"); return 0; } sum = 0; for (j=0;j<i;j++) sum += Up[i][j]*Y[j]; Y[i] = (VpB[i] - sum)/Up[i][i]; } /* Triangle Superior */ if (!U[4][4]) { printf ("Division by zero !\n"); return 0; } X[4] = Y[4]/U[4][4]; for (i=3;i>=0;i--) { if (!U[i][i]) { printf ("Division by zero !\n"); return 0; } sum = 0; for (j=i+1;j<5;j++) sum += U[i][j]*X[j]; X[i] = (Y[i] - sum)/U[i][i]; } printf ("%.2f + %.2fx + %.2fx2 + %.2fx3 + %.2fx4\n", X[0], X[1], X[2], X[3], X[4]); return 1; }

This will give a 5th order polynom which is fairly enough to describe the graphs we are dealing with. Just replace the precomputed curve by the values found in the table :-)

Boomer - 14-5-2006 at 09:12

Hi both!

You asked for a formula to calculate det pressure (CJ-pressure). The simplest is

P = d*D^2 / 4

with P in gigapascal, d in g/cm^3 and D in km/s. It gives very good values, e.g. 33 GP for RDX at D=1.77 (measured 34GP), and around 39GP for HMX. Using m/s and .../400 gives P in kilobar.

It is derived from

P = d*D^2 / (y+1)

with y = a factor derived from the latent heat of the reaction gasses (I will look it up more precisely). y is close to 3 for most CHNO explosives, hence the .../4 in the above formula. (Both from Cooper)

I also have an old formula at home, using density of HE, density and streaming velocity of reaction products etc. It is more complicated and gives values 20% too low. But it is the one used in old books, where TNT was said to have only 160 kilobar. Again, I will look it up. This old one is from one of the Stettbacher books IIRC.

Good work btw!

[Edited on 16-5-2006 by Boomer]

sylla - 15-5-2006 at 09:12

d in g/cm^2 ? I guess you mean cm^3...

enhzflep - 17-5-2006 at 23:27

RE: pressure formula.
Great, thank-you very much. Haven't checked to see if it's the same one
currently implemented, but if not I'll add it in along with appropriate
selection control. Don't worry, I knew what you meant ;)

RE: code
Any chance of a re-post of that last code example, as an attachment?
The board's formatting software has screwed up the code and removed
serveral sections. Some of which i possess insufficient intellect to correct.
Ignore U2U Regarding same.

Additionally, have nearly completed a brand spankin' new version,
have already implemented suggestions above including -

*text file output
*bracket implementation
*configurable units i.e CC/CI - g/oz - KJ/KC
*library of all entered explosives
*ability to oxygen balance (no CO/H2O yet)
*ability to do VDet of mixtures.
*use of 1 font across program to avoid invisible parts of boxes.

Have only had trouble implementing different colour scheme and ability
to size the dialog boxes. This entails quite a bit of really heavy duty
programming, I'm prepared to do that but cannot find refferences on how
to go about it at this stage.

I've decided against using the system registry for configuration data
instead opting to write a ".cfg" into the operating directory. This will mean
that the program is much easier to remove entirely from one's system;)

Now, just gotta write that damned help-file :( What a pain in the ass.

Any further feedback, further requests or comments entirely welcome

sylla - 21-5-2006 at 09:04

Here is the source code. There are a few odds with the units selection but it should work fine for testing

Attachment: ChemCalc7.rar (10kB)
This file has been downloaded 673 times