Panair
Polishing FDM with Panair, Blender, and LibreOffice Calc
Best of Flight Gear models at now are going with JSBSim FDM, based on real data from pilots manuals. It's most conventional way to make realistic model. But sometimes manual is not accessible - then commonly Datcom is used to advance initial Aeromatic data to appropriate state.
Digital Datcom is USAF stability program from 60x, based at theoretical formulas and polynomial approximations. It not deals with real model, but with parameters as perimeter or diameter. Anyway, for most of cases that's enough - I compared its results with MiG-15 manual data at time, and did not found any difference.
But sometimes it's not so. For example, complex fuselage of Su-15 could be calculated with Datcom way not better than one demonstrated at screen. Also, You may notice that form of tail differs. That would be fixed, but by two parts tail only. Datcom principally can't calculate wing with smooth swept angle changing.
Some would try to use CFD then, as NASA OpenVSP. Especially ones having brand new videocards, since at now its VSPAero processor testing CUDA capabilities. Some, as me, don't have such capabilities - at my ten years old PC calculation of one point, of needed hundreds, took hours at example F-16 model with it, prior I shut it down of course.
For tests I chose Boeing PanAir program instead. That tool from late 70x do have some attractive advantages:
- It's completely freeware, with Fortran code opened. So, theoretically, if something goes wrong, You could look into it and fix something.
- It has ability to calculate supersonic cases, which other implementations of similar approach, as modern and opened XFLR5, don't.
- It's from 70x and 80x, uses panel method, means that instead of full 3D CFD task it calculates flat task really; so should work fast.
Two last assumptions turned out to be true more or less, so, let's get trough FDM calculation process step by step, to figure out how use it; PanAir itself meant to be downloaded and compiled by gfortran already.
Step first: FDM 3D model preparation
You need to make model, as shown at screens. Example Su-15 model is presented as .ac file at downloadable pack by link below, with all other needed additional stuff.
It's got to be:
- ~500 polygons approximately model. Presented code seems to accumulate some error with calculations - twice more detailed model takes twice greater time and produces twice wronger results.
- Half of closed body. Actually, full body model can be calculated too, but it's just unneeded waste of time.
- Made pixel to pixel, without any holes or slides. PanAir, as said, allows to calculate models where, for example, points of wing root are at line of gap at fuselage surfaces, differs than points of fuselage. But at most of cases attempt to do so produces thing which would be called, by analogy, "abutment hell" - You would not want to know what it is.
- With sizes at feets. Model with sizes at meters gives wrong results.
- With more or less even grid.
- With normals changing less than 20deg at each part, and each side of two-sided surfaces as wing or elevator.
- Without points at negative areas of "X" axis.
- Which could be divided at parts of pure rectangular form in means of points quantities, or started as cone with rectangle at just next point, or finished as cone, or both, but not of any other forms.
Otherwise program just'll quit with an error.
Step second: Cut
Cut it to parts with short simple names, again, as shown at screens. Note that wing and elevator are cut to tip and one torus surface.
Step third: Tables Export
For export I made GPL Blender plugin at Python, code of that is presented at downloadable pack. You need to install that prior, of course - as You, I suppose, installed .ac files export/import plugin earlier.
It's not a product, but working tool - without safety shields for most of possible errors cases. So I would recommend to launch Blender from console, where You would see output log, and could stop it by Ctrl-C at case of hanging.
Plugin has some options, which should be set correctly prior exporting. "First Axis" is sequence at which points are going at rows of resulted table of export - it's got to be "X" for all surfaces except main surfaces of wing and elevator. "Second Axis" is sequence at which points are going at columns - it's got to be "Z" for most.
Exceptions are: main surfaces of wing and elevator, and tip of tail.
For first two "First Axis" got to be "-X" and "Second Axis" should be equal "Y"; also, "Two-sided wing" option gotta be on. Main surfaces are declared at PanAir as toruses, starting from back, going clockwise at rows, and from fuselage to tip at columns. Note that first point got to be found correctly - if wing has trapeze form than "First point" parameter could stay "First Axis", but for swept wing You got to set "Second Axis", as points of tip are more distant from nose by X axis than points of root.
Axises for tip of tail are "X" for first and "Y" for second.
It should be mentioned also, that plugin can make some other peculiarities, caused by not so straight input data. For example, first point of back of fuselage at presented example model is not lower one, so attempt to export it correct way, "X/Z", leads to an error, while with "X" and "-Z" it's exported correctly, and PanAir gets that inversed order without errors.
Step fourth: Input File Combining
Here goes some code. PanAir input file commonly consists of next sections, which are explained at base of presented Su-15 model example.
$title simple Su-15 full configuration test
Seems to be got to have two strings at least - does not matter much what's going on here.
$datacheck 0.0
0.0 - main program run, 1.0 - full data check, 2.0 - short data check. Should be 0.0.
$symmetry xz plane =misymm mjsymm 1. 0.
That configuration is for half of body model. For rotated surface tests, which are explained later, got to be 0. 0. .
$mach number =amach 0.3
Self explanatory.
$cases 1.0
Actually program could calculate up to four cases with close enough, not more than ~5deg from middle, alphas. That option is not used, so 1.0.
$angles of attack =alpc 10.00 =alpha 10.00
For one case alpc=alpha.
$references for accumulated forces and moments =xref yref zref 40.9 0. 0.22 =sref bref cref dref 212.74 28.08 12.92 95.0
xref, yref, zref are coordinates of center of mass. For example, at Su-15 these was not known, but fact that CM got to be forward of main wheels, and as close to aerodynamic center of wing, which is ~0.33 of mean aerodynamic chord for triangle wing at subsonics, fact which is known from literature, limited possible range to that data. With other CM it would just sit at tail, or could not rise nose up at takeoff.
sref is surface of whole wing at XY plane, both left and right sides, without part which is inside fuselage. bref is full span of wing, both sides. dref is length of whole model.
$points of body =kn cpnorm 3. 2.0 =kt 1. =nm nn netname 5. 6. ascanopy =x(1,1) y(1,1) z(1,1) x(*,*) y(*,*) z(*,*) 15.4378 -0.0000 2.8555 18.0107 0.9549 3.0279 ...
kn is number of meshes at section. cpnorm - 2.0, or 0.0, is when normals of panels, means polygons of mesh, are taken before automatic surface to surface fit attempt. 1.0 - after, but, as initial model is already fitted point to point, that does not matter. 2.0 should be calculated faster. kt - boundary condition. You would want to know what it is from documentation. For subsonic body and wing it's 1., for back hole of engine or body 5., and so called wakes are explained later.
=nm...
- here goes Your exported data. Note that "=" is commentary sign, whole string with it is omitted from calculation.
I would recommend to export by different sections not all and back, but fuselage, wing, tail, elevator - as at example. Then You would check results, which are produced by different parts then, at panair.out file separately.
$trailing wakes from wing and elevator =kn cpnorm 2. =kt 18. =inat insd xwake twake netname aswing 1. 90. .0 wingwk aselev 1. 90. .0 elewk
$trailing wakes from body =kn cpnorm 2. =kt matcw 18. 1. =inat insd xwake twake netname asfdown 3. 90. .0 blowk asfup 3. 90. .0 bhighwk
That part could be a bit tricky. Way PanAir calculates model demands wakes to be defined - without these output data are worser. Wakes starts from back side of surfaces and fuselage at float, and goes to ~1.5 of fuselage length from nose. insd parameter defines which side of surface is that end actually; trickness is that defines in reversed order.
So, if root of wing is really first string of data, last is tip, and at string it starts from greater Xs, then it means 1-root, 2-forward side, 3-tip, 4-back side; but then this sequence is inverted, so 1-back side, 2-tip, 3-front, 4-root.
For tail and body initial sequence is 1-low, 2-back, 3-up, 4-front, and used inversed one is 1-front, 2-up, 3-back, 4-down. Dude who made it, You know, was programmer. Anyway, it's how that is at official example and docs.
Note that tail wake is absented - for some reasons program can't calculate it, while there is not time to figure out why. It works as it works, somehow, and this is good already - it would not work completely.
$eat - liberalized abutments =epsgeo igeois igeout nwxref triint iabsum 0.025 0. 0. 0. 0. 1.
It's data for automatic abutment processor. Program has other one, by which You could say to PanAir what connected with what manually, but remember inversed order and forget about it. epsgeo parameter defines how far collided points could be in means of smaller chord of polygon of model, and it shoul't be too big. iabsum defines what part of abutment processor logs should be printed at resulted data, default value is ok.
$forces and moments summary for nonwake networks
That's what we need.
$trefftz plane analysis
By experience, original method produces bad data sometimes, especially at oversonics and greater alphas. For these cases optional Trefftz analysis is better for sure, actually it should be default one, but let's talk about it a bit later.
$printout options =isings igeomp isingp icontp ibconp iedgep 0. 0. 0. 0. 0. 0. =ipraic nexdgn ioutpr ifmcpr 0. 0. 0. 1.
ifmcpr=1. means we get not only sum forces at output file, but for each surface.
$end
End means end.
Step fifth: Calculation
Of course, You can run PanAir manually, typing input file name each time, and deconstruct output file manually too. Yes, You can. But better way is to use proposed script, at Python, since Python is already installed with Blender. Call string is
python alpha-mach-coeff-test.py input_file.inp
meant painar exec file is at same directory with input file. Results are at input_file.out, also, these are printed at console. If there is no results, then You made something wrong, look panair.out file. Testing Mach numbers and attack angles are declared at beginning of script - just change these arrays if You need other ones. One point calculation takes around twenty seconds for example model at ten years old PC, so whole case is question of minutes commonly.
Step sixth: Data Mining
That part is most interesting. I would not recommend You to use R right now, not yet, but You would need some LibreOffice Calc most probably at least. Calculated data is at input_file.out - import it at Calc.
What You may see from beginning - resulted cl, lift coefficient, and cdi, drag, are completely wrong for presented Su-15 model. Values as -37.902 for drag, You know. So, moments are wrong completely too.
But - bright side - results of additional Trefftz method, tre_cl, tre_cdi, these are more or less nice - see comparison with Datcom results at presented graphs for Su-15 at 0.3Mach case.
Clearly it does not care about critical angle and so - still, with complex bodies, as "Su-37" or "Su-50", and other wingbody stuff, thing as that could be only option.
Also, I would, add that data for subsonics, ~1.0Mach, are wrong with it, great angles of attack are omitted at higher Machs, so linear
FORECAST()
function of Calc is Your good friend - learn it to fit wrong or absented data. I tried to figure out how make PanAir calculate 3.0Mach, as it proposed working up to four, but seems to be some additional things needed here, look "Boundary Conditions" of official manual for that. I failed due to absence of time, maybe You'll make it - if it really can be made, of course. But if You are at same situation as me then these higher Machs values could be approximated too. Luckily, Su-15 is limited by 2.0Mach anyway.
Step seventh: Moments
As that thing failed at normal moments coefficients calculations, You need to recalculate it from lift force coefficient, by known fact that
cm=cl*shoulder_from_CM_to_CAERO/cref
At "Su-15" CAERO is known as 1/3 of mean aerodynamic chord, while CM was received by tests; result is shown at screen.
Step eight: Side forces and moments
Screen is self explanatory.
Step nine: Controlled Surfaces
Example and result are at screens.
Theoretically, PanAir could calculate effects at controlled surface as part, oversonic shadowing, so on, but at fact results are quite bad or absented at all of these cases, and only option is to calculate rotated surface alone, without anything else. Code for rotation is
$rea of elevator roation =nrearr 1. =ntr 1. =k1 k2 1. 3. =x1 y1 z3 x2 y2 z2 60.1421 3.6781 0.6448 60.1421 4.6781 0.6448 =phi -8.00
where k1..k2 are numbers of meshes at list, from which to which all is rotated. Then goes coordinates of two points of rotation axle, and then angle. That part got to go after elevator points declaration of course, and prior wake. Calculating script and example are attached; note that at such simple case main cl coefficient is more or less correct, and could be used, but cm better to be recalculated by, as simple look at aero file tells
cm=cl*CM_to_elevator_axle_shoulder/bref
Since not alpha changed, but phi, other script is used:
python rotated-coeff-test.py input_file.inp
but principle is same.
Note that for simple sole surface main program method works, and produces more realistic data than Trefftz.
Resume
So, what we got there? Look at graphs carefully. Tool which was made and used for years by leader company at market, with, at time, supercomputers, which, as NASA says at papers, calculated 3.0Mach with 1% of error, produces quite less exact results than outdated military stuff, at quite common 0.3Mach case. That can't be, so, most probably, something was lost in translation.
Also, I could add, program compiled with all debug drops out with some array errors. So, maybe, matter is modern Fortran differs from one used initially at some means. If You think of own self as of smart, You could try to fix something at 60000 strings of code. If, let me say it clearly, they did't cut something out prior making program public, You would even achieve supertool, which current version isn't.
At my view, to test some other approaches instead would be correct idea; anyway, that thing could be used at some cases - first case which goes to mind is complex body form at lower Mach/attack angles, but greater than some glider has. Also, to just export simplified model with minimal changes, to run some script and get resulted table - in minutes - could be more attractive for some developer than manually calculate all needed data for Datcom. Even by price of some simplicity of output.
And - last but not least - it would be quite useful if You doubt in Datcom results and would want to polish it a bit, as I made with Su-15 model.
External links
Example PanAir cases made by NASA, note, for older version of program
Su-15 Panair Pack, including Blender export script, run Python scripts, .ac and .inp examples
Author
Victor Slavutinsky, vitosnet at mail dot ru