C++ Question organizing code

[malisle]

I am currently in the process of writing some functions and code snippets in C++ for an yet unspecified add-on and I have some questions about the ways to organize the code. I made a few functions that serve as first approximation of a solar powered probe power system (shadow cones, battery charge/discharge cycles,solar panel efficiency as a function of temperature etc.) and basic thermodynamic system (hull and panel equilibrium temperatures), and I find myself constantly recalling the same few functions. For example, function that checks if the vessel is in the shadow is used both for the thermodynamic and power system, and is called minimum two times each step. When I executed it directly in opcPreStep() and assigned a return value to some placeholder variable visible to other functions number of executions per step was reduced to one, whether it was needed or not. So, my question is how do you keep the amount of unneeded code executed in each step to a minimum? What is the best way to organize interaction between "lower" and "higher" functions? Any other advice on the ways to organize and optimize the code would be appreciated.

 

[Hlynkacg]

Can't speak for others. I'm a bit of a c++ newbie my self, but for me what works best is to create lower order functions that return a value and then call that value in your higher functions.

To use your example.

For instance the following function determines if we are in sunlight...

// Figure out if we are in sunlight
bool LM_Enviroment::InShadow()
{
	VECTOR3 rBodyPos, gsunPos;	// Vectors to be used.

	// Calculate Shadow cone
	oapiGetGlobalPos (rBody, &rBodyPos);								// get global position of reference body
	double bodySize = oapiGetSize (rBody);								// get Radius of reference body
	double rBody_Sun = dist (_V(0,0,0), rBodyPos);						// get distance of reference body from sun
	long double vessel_Sun = dist (_V(0,0,0), gPos);					// our distance from the sun
	long double Lcone = (rBody_Sun * bodySize) / (6.96e8 - bodySize);	// Shadow cone = (distance from sun * rBody radius) / (Sun radius - rBody radius)
	long double diff	= vessel_Sun - rBody_Sun;						// Difference between rBody's distance from the Sun and ours (is it between us?)
	long double radius	= length (rPos);								// Orbital radius (our distance from rBody's center)
	long double angle	= atan (bodySize / Lcone);						// Angle of shadow cast by rBody
	long double alpha	= pow (radius, 2) - pow(diff, 2);				// 
	long double beta	= sqrt (alpha);									//
	long double gamma	= Lcone - diff;									//
	long double angle_V	= atan (beta / gamma);							// 

	// Return a value
	if ((angle_V < angle) && (diff > 0)) return true;
	else return false;
}

once you declared it. Then in your higher functions can call it using something along the lines of...

if (InShadow) do X;
else do Y;

PS if you've figured out temperature equilibrium drop me a line, my Temperature code keeps giving me wonky results.

 

[ADSWNJ]

I would commend you for thinking about performance and efficiency at this early stage of writing your code. So many times, this is a topic ignored until the users complain about performance, and even then the typical developer response is more RAM or faster CPUs needed .

You should think about the overheads of function calls, and particularly the memory allocations and parameter handling for each one. Use consts whenever you don't need to update things. Think about making the definition inline (trading more code for less calls). Finally, see if you can figure out where the code is spending the bulk of it's time, and focus your efforts on optimizing calls and loops in there. (For Pro coding, there's lots of products that assist in this space.)

For orbiter, just think about what update loops and where you need to call for a recalc. For example - if your MFD update interval is 0.05 secs, and you have a heavy calf, then do you need to redo it every step, or could you optimize it with your own sim timer trigger?

 

[Enjo]

It's best to organize your code into many hierarchies, two is the minimum, where Hlynkacg's example would be the lowest level. If the lowest level can be made using static methods, it's already good. At some point the higher level will need to store state in variables - make classes for that.
Look into Launch MFD's sources to get an overview how it could be organized. In the lib/Math folder you'll find many examples of the lowest level, environment independent classes, which contain more or less static methods, but also classes which are variable dependent, like SpaceMathBody class, which requires body's mass and radius for further calculations, and stores them internally. The lib/Systems/Vect3 class contains more of a "fat-interface" example, where the vector's length is a method of the Vect3 and operates on its data.

Check how all these low level methods are called in Launch MFD, using some IDE (look for references of methods), and try to measure how many lines of code is saved, using this simple technique.

Now to refer to ADSWNJ's reply, the Vect3::len() method should have been declared inline. This wouldn't hurt for sure, because its pretty small, but don't make any guesses for other non-obvious cases. At some point you have to run your code through the profiler, before making any moves, that could be a simple waste of time. And never sacrifice code readability for potential speed up. I'm with the "Buy a better CPU" guys here, unless you've profiled the code and it turned out to be an obvious mistake, like using a vector in a place where you should have been using a map, whose search time is logarithmic vs. linear in case of vector.
As the #C++ channel's saying goes: "Make things work, profile, optimize. In that order"

Because it's hard to profile an MFD for Orbiter, I've always thought that my autopilots take a lot of CPU time, because they contain a lot of math. Once I've started using DX9, it turned out that the inline graphical client's video refresh is the bottleneck of Orbiter, not my APs. Math stuff is calculated pretty fast on today's CPUs, especially on 64-bit OSes and applications. You can bet on a 2.5 speed up for floating point operations there (not the case for Orbiter though), for something that's already calculated pretty damn fast. It's the I/O operations in general (graphics, file) that take most of the CPU time, and some obvious mistakes, like in the example with map and vector.

For code organization and general C++ coding, I recommend the article here:
http://www.orbiterwiki.org/wiki/Developing_modules_in_C++_(basic)

 

[MeDiCS]

Programmers waste enormous amounts of time thinking about, or worrying about, the speed of noncritical parts of their programs, and these attempts at efficiency actually have a strong negative impact when debugging and maintenance are considered. We should forget about small efficiencies, say about 97% of the time: premature optimization is the root of all evil. Yet we should not pass up our opportunities in that critical 3%

I think this is one of the few 'rules' a serious programmer should take to heart. The question is whether you know which part of your code fit into those 3%. (hint: timer and profiler)

With that said, there are optimizations that actually do make code more readable and simpler. You shadow test is one such example. Lets assume your code is something like:

void UpdateThermo(...)
{
    ...
    if (InShadow())
    ...
}

void UpdatePower(...)
{
    ...
    if (InShadow())
    ...
}

DLLCLBK void opcPreStep(double simt, double simdt, double mjd)
{
    ...
    UpdateThermo(...);
    UpdatePower(...);
    ...
}

We have two calls to InShadow() in the same step, and they will, as a rule, return the same value. Now, to realize that they return the same value takes a bit of thinking, and that (relatively important fact) may be missed by some reading your code. If you do cache the return value instead:

void UpdateThermo(..., bool shadowed)
{
    ...
    if (shadowed)
    ...
}

void UpdatePower(..., bool shadowed)
{
    ...
    if (shadowed)
    ...
}

DLLCLBK void opcPreStep(double simt, double simdt, double mjd)
{
    ...
    bool shadowed = InShadow();
    UpdateThermo(..., shadowed);
    UpdatePower(..., shadowed);
    ...
}

You not only reduce redundancy, but also give a clue to the reader (which may be yourself) that the vessel will either be in the shadow or not, and that that will not change until the next step. So, my advice is that you should keep code as simple and with as little redudancy as possible, while remaining sane, and reserve 'optimization efforts' to the few parts of your code that actually needs it.

 

[malisle]

Thanks for all the answers! Some of the logical answers I can see by myself but solutions that require a more detailed knowledge of C are still somewhat foggy. Something to keep me entertained for a while

PS if you've figured out temperature equilibrium drop me a line, my Temperature code keeps giving me wonky results.

Here they are. Second one checks if the vessel is in shadow, first one is needed by third, and the third calculates the solar panel temperature. Temperatures will be as expected, but their dynamics seems to be linear instead of aperiodic when increasing/decreasing due to day/night cycles.

double shipname::f_irradiance (void) [COLOR="SeaGreen"]//returns [W/m^2] of sun radiation at the current position of the ship[/COLOR]
{
	VECTOR3 vesselpos;
	GetGlobalPos(vesselpos);
	return ((pow(6.96e8,2)) / (pow(length(vesselpos),2)) * 73e6);
}

void shipname::f_insunlight (void)  [COLOR="SeaGreen"]//checks if the vector that connects ship and sun passes through reference body[/COLOR]
{
	OBJHANDLE sun=oapiGetObjectByName("Sun");
	OBJHANDLE mainbody=GetGravityRef();

	VECTOR3 vesselpos, mainbodypos, sunpos,v1,v2;
	double mainbodyR=oapiGetSize(mainbody);
	double d;
	
	oapiGetGlobalPos(sun, &sunpos);
	oapiGetGlobalPos(mainbody, &mainbodypos);
	GetGlobalPos(vesselpos);

	if(mainbody==sun){insun=1;} 
	else if (mainbody != sun && length(vesselpos) < length(mainbodypos)){insun=1;}
	else if (mainbody != sun && length(vesselpos) > length(mainbodypos))
	{
		v1= _V(sunpos.x-vesselpos.x,sunpos.y-vesselpos.y,sunpos.z-vesselpos.z);
		v2= _V(vesselpos.x-mainbodypos.x,vesselpos.y-mainbodypos.y,vesselpos.z-mainbodypos.z);
		d=length(crossp(v1,v2))/length(v1);
		if(d>mainbodyR) insun=1; else insun=0;
	}
}

void shipname::f_solpan_temperature (double simdt) [COLOR="SeaGreen"]//solar panel temperature [°C][/COLOR]
{
	double irradiance = f_irradiance();
	double sunangle   = (1-abs(orient.y))*(1-abs(orient.x));

	double eqilibium_temp = pow(((SOLAR_PANEL_ABPSORTIVITY * SOLAR_PANEL_AREA * irradiance * insun) / ((SOLAR_PANEL_EMISSIVITY_FRONT * SOLAR_PANEL_AREA + SOLAR_PANEL_EMISSIVITY_BACK * SOLAR_PANEL_AREA) * STEPHAN_BOLTZMAN)),0.25);
	static double temperature = eqilibium_temp;

	double absorbed_heat = irradiance * SOLAR_PANEL_ABPSORTIVITY * SOLAR_PANEL_AREA * sunangle * insun;
	double emitted_heat =  (SOLAR_PANEL_EMISSIVITY_FRONT * SOLAR_PANEL_AREA + SOLAR_PANEL_EMISSIVITY_BACK * SOLAR_PANEL_AREA) * STEPHAN_BOLTZMAN * pow(temperature,4);
	double dT_dt = (absorbed_heat - emitted_heat) / (SOLAR_PANEL_MASS * SOLAR_PANEL_SPECIFIC_HEAT_CAPPACITY);[COLOR="SeaGreen"]//from Q=m*c*dT[/COLOR]
	temperature+=(dT_dt*simdt);
	[COLOR="SeaGreen"]//sprintf(oapiDebugString(),"irradiance[W/m^2]=%.2f, eqilibrium T[K]=%.2f, absorbed Q[W/m^2]=%.2f temperature[°C]=%.2f",irradiance,eqilibium_temp,absorbed_heat,temperature-273.15);[/COLOR]
}

 

[Urwumpe]

One small optimization for you there:

static const double REF_DISTANCE = pow(6.96e8,2)/73.0e6;
GlobalPos(vesselpos);
return REF_DISTANCE / pow(length(vesselpos), 2)

While this does exactly the same as your code, this one can be easier optimized by the compiler. In your version, you depend heavily on the performance of the optimizer to notice that you have a constant part and a variable part in the term. By explicitly splitting it into constant and variable part, the compiler can find a better solution.

Next optimization:

pow(length(vesselpos), 2) is equivalent to dotp(vesselpos, vesselpos).

Saves you two of the worst floating point instructions and makes your program execute faster... not that it really matters... but it can and it is a good easy example how you can help the compiler to help you.

Of course: Still write proper comments around the optimized code, so you can understand it, ideally write them to complete, that you could understand your code after 2 years of not looking at it and while being completely drunk. In 2 months already, such good comments mean that you need only 10 seconds of reading, where you would without comments need 5 minutes to find out what the hell this function does and how REF_DISTANCE is defined.

 

[martins]

If you know up-front that the shadow check will be performed at least once per frame, doing the calculation in opcPreStep and storing it for later use, as suggested above, is a good idea.

If however it is possible that the check is not performed at all in some frames, this would still be wasteful. In that case, it would be better to use a flag to indicate the validity of the shadow check, reset it at the beginning of each frame, and recalculate the shadow check only on demand, if the flag indicates that it is not valid. Something like this

void MyVessel::clbkPreStep (...)
{
   shadow_check_uptodate = false;
}

bool MyVessel::InShadow()
{
  if (!shadow_check_uptodate) {
    in_shadow = ShadowCheck();
    shadow_check_uptodate = true;
  }
  return in_shadow;
}

The (potentially expensive) ShadowCheck method is then called at most once per frame. ShadowCheck() should be private, so that there is no chance somebody calls it unnecessarily from outside. (as of course should be shadow_check_uptodate and in_shadow).

 

[csanders]

Here's a method giving the shadow state its own object. This is similar logic to how I commonly cache images and such as.

//ShadowState.h
class MyVessel; //This is so the compiler won't balk at the MyVessel * object.
class ShadowState
{
public:
    ShadowState(MyVessel *vesselPtr);
    bool InShadow();
    void Reset();

private:
    bool upToDate;
    bool inShadow;
    MyVessel *vessel;
    bool ShadowCheck();    
};

//ShadowState.cpp
#include "ShadowState.h"
#include "MyVessel.h." //this is so the code can actually use the MyVessel *

ShadowState::ShadowState(MyVessel *vesselPtr):
   upToDate(false), inShadow(false), vessel(vesselPtr)
{
}

bool ShadowState::InShadow()
{
    if(upToDate)
        return inShadow;

    inShadow = ShadowCheck();
    upToDate = true;
    return inShadow;
}

void ShadowState::Reset()
{
    upToDate = false;
}

bool ShadowState::ShadowCheck()
{
    //some code that checks for shadow.
    //presumably uses the MyVessel *

    //THIS IS THE 'EXPENSIVE' FUNCTION TO CALL
}

//MyVessel has an ShadowState object called 'shadowState'
//FYI: MyVessel's constructor will look something like this:
//MyVessel::MyVessel(...):shadowState(this){}

void MyVessel::clbkPreStep (...)
{
   shadowState.Reset();
}

//Now these functions won't call the expensive 'ShadowCheck()' function
//unless it needs to be called. And it won't be called more than once between
//resets.
void MyVessel::UpdateThermo(...)
{
    ...
    if (shadowState.InShadow())
    ...
}

void MyVessel::UpdatePower(...)
{
    ...
    if (shadowState.InShadow())
    ...
}

 

[martins]

Incidently, this could be pushed even further. If you are defining the shadow condition as a binary state, I assume that you are not fussed about partial illumination at the light-shadow transition. Since this transition is not instantaneous, due to finite size of solar disk and atmospheric effects (if applicable), checking at each frame is probably overkill. You could therefore reduce the check frequency:

void MyVessel::clbkPreStep(simt)
{
  if (fabs (simt-last_check_t) > update_interval) {
    shadow_check_uptodate = false;
    last_check_t = simt;
  }
}

If you set update_interval to 1s, then at 100fps you would skip 99% of the shadow update calculations.

 

[Enjo]

unless you've profiled the code and it turned out to be an obvious mistake, like using a vector in a place where you should have been using a map, whose search time is logarithmic vs. linear in case of vector.

There's a bit more to say in this topic actually, and it's really the place where discussions about speed make much sense, because for larger amount of computations it can make a difference of runtime of one hour vs. one year or more.

To show an example, last time we were discussing with ADSWNJ the following problem: Imagine that you're in a circular orbit and you want to deorbit your craft. You're interested in knowing where your potential landing site would be if you applied, say a -200 m/s retrograde burn at current frame. Therefore, you need to learn the 3D position vector at the place, where your orbit's position vector's length is equal to the radius of planet (in other words your orbit crosses the surface). You'd then convert this position vector to latitude and longitude using Orbiter API and viola.

Now how to get this position. My proposal was to use [ame="http://orbithangar.com/searchid.php?ID=3825"]KOST[/ame]. Looking at it from a broader perspective, KOST is a configurable (given initial state vectors) function, that translates time into position, so pos = KOST(t), where time is the amount of time passed from the moment of providing state vectors. As a side note, the state vectors' velocity component needs to be shortened by 200 m/s, as already said, before performing the KOST simulation. And here we have a typical problem. You don't know the optimal time which provides the required pos, ie. such that it crosses the planet's surface. You only know that the maximal time range is in (0, orbitPeriod/2), because the periapsis is at orbitPeriod/2, and if you don't reach surface then (|pos| > Radius), then it means that you'll never reach it with the dv of -200 m/s, and you need to decelerate more. But let's assume now that you do reach the surface in this time range.

Amateur's approach (Naive algorithm):
Iterate KOST function with a constant time step, like 10s and check on each iteration if (|KOST(t_i)| < Radius). If it does, you return pos = KOST(t_i). If the orbit's period was 6000, then T/2 = 3000, so pessimistically you needed 3000 / 10 = 300 iterations to solve this simple task. Let's assume that it turned out that 10s precision was giving you unstable results, and you need to increase the precision to 0.1s. You get 3000 / 0.1 = 30000 iterations! You'd surely start feeling this in your simulation. This was an example of linear complexity, or O.

Pro's approach (Divide and Conquer):
Since you know the maximal time bound, which is T/2, you start off from calculating the KOST(t) for T/2 and note the sign of (KOST(t) - Radius), it should be negative in this case, because we said, that -200 m/s is enough. If it isn't, break the calculations, because you'll not reach the surface. Now slice the current reference time, T/2 into 2 partitions, and search for KOST(t) in the middle of the two partitions, so for t = T/2 - T/4 = T/4 and note the sign. If it has changed, then it means that the solution is in the furthest partition (between T/4 and T/2), and if it hasn't, then it's in the first partition (between 0 and T/4), let's assume that it's in the 2nd partition (realistic case). We slice it into next 2 partitions, and proceed into the middle of them, so t = T/4 + T/8 = 3T/8 and test KOST(t) there, and so on, until we reach the stated precision (difference between subsequent time slice values). The total number of iterations will be 2*log(maxArg / precision), so in our case 2*log(3000 / 0.1) = 21! This will be barely noticeable in your sim, because it's just a little bit of pure math, and no I/O operations. The general idea is to slice the problem into many binary sub problems, where the answer to each one is either true or false. This is an example of a binary search with logarithmic complexity, or O(log), and this is where you should be spending most of your optimization effort, when in your app you need to guess an optimal argument of function numerically. You could do it analytically, but in general that works only in schools on unrealistically simplified examples. Naturally, if you're able to solve such complex problems analytically, then hats off.

I've used this technique in Direct Ascent recently, where my function is PEG, from which I try to acquire the engine thrust, that gives me a concrete time of orbital insertion, that is such that is equal to target's arrival on my insertion position, so deltaArrivalTime = PEG(engF). The delta time must be 0 in my case, not Radius, as in previous case, but it only means that in my case the reference value is 0, while in the previous case the reference value was Radius.

This problem is so generic, that I've decided to write a generic binary solver for such problems:
BinSearchArg.cpp
BinSearchArgSubject.hpp