Search code examples
copenclssesimddifferential-equations

How to optimize valve simulation logic?

This is an simple logic-programming and optimisation etude, that I've created for myself, and somewhat stumbled upon it.

I have a numerical simulation of simple scheme. Consider some reservoir (or a capacitor) Cm, which is constantly pumping up with pressure. Lets call its current state Vm:

enter image description here

At its output it has a valve, or a gate, G, that could be open or closed, according to following logic:

  • Gate opens, when pressure(or voltage) Vm exceeds some threshold, call it Vopen: Vm > Vopen
  • Gate remains open, while current outrush Ia is greater than some Ihold: Ia > IHold
  • Gate conducts power only out of reservoir (like a diode)

I am doing numerical ODE solving of this, i.e. determining Vm and Ia at each (equal, small) timestep dt. Have three variants of this:

Variable types:

float Vm=0.0, Ia=0.0, Vopen, Ihold, Ra, dt;
int G=0;

Loop body v1 (serial):

Vm = Vm + (-Ia*G)*dt;
G |= (Vm > Vopen);
Ia = Ia + (Vm*Ra*G)*dt;
G &= (Ia > Ihold);

Loop body v2 (serial, with temp var, ternary conditionals):

int Gv; // temporary var
Vm = Vm + (-Ia*G)*dt;
Gv = (Vm > Vopen) ? 1 : G; 
Ia = Ia + (Vm*Ra*Gv)*dt;
G = (Ia > Ihold) ? Gv : 0;

Loop body v3 (parallel, with cache):

// cache new state first
float Vm1 = Vm + (-Ia*G)*dt;
float Ia1 = Ia + (Vm*Ra*G)*dt;
// update memory
G = ( Vm1 > Vopen ) || (( Ia1 > Ihold ) && G);
Vm = Vm1;
Ia = Ia1; // not nesesary to cache, done for readibility

the G was taken from building up the following truth table, plus imagination: enter image description here

Q:

  • Which is correct? Are they?
  • How third variant(parallel logic) differs from first two(serial logic)?
  • Are there more effective ways of doing this logic?

PS. I am trying to optimize it for SSE, and then (separately) for OpenCL (if that gives optimisation hints)

PPS. For those who is curious, here is my working simulator, involving this gate (html/js)

Solution

    1. By the overall description they are the same and should all fulfil your needs.

    2. Your serial code will produce half-steps. That means if you break it down to the discrete description V(t) can be described as being 1/2 dt in front of I(t). The fist one will keep G changing in every half-step and the second one will synchronise it to I. But since you are not evaluating G in between in doesn't really matter. There is also no real problem with V and I being a half step apart but you should keep it in mind but maybe you should use for plotting/evaluation the vector {V(t),(I(t-1)+I(t))/2,G(t)}.

      The parallel code will keep them all in the same time step.

    3. For your pure linear problem the direct integration is a good solution. Higher order ode solver won't buy you anything. State-space representation for pure linear systems only write the same direct integration in a different way.

    For SIMD optimisation there isn't munch to do. You need to evaluate step-wise since you are updating I by V and V by I. That means you can't run the steps in parallel which makes many interesting optimisation not possible.