/* B-2; version 1.0 */
/* (c) copyright Ian Hickson 2000, distributed under the GNU GPL */

#include "dynamic-2d-arrays.h"
#include "lists.h"
#include "runge-kutta.h"

/* Circuit Constants */
#define C 4.7e-6
#define R 1e3
#define L1 20e-3
#define L2 10e-3

/* Details of range */
#define Start 0
#define End 20e-3
#define StepSize 1e-5 /* MUST be suitably small as runge kutta is not stable! */

/* Differential Equation: 

   d3i    1   d2i   1    1     1   di     1
   --- + ---- --- + -- (--- + ---) -- = ----- e(t)
   dt2   L1RC dt2   L1  L1C   L2C  dt   L1L2C
         \_b/       \______c_____/      \____d___/
 */

nlNumber e(nlNumber x) {
  return x < 5e-3 ? 1.0 : 0.0; 
}

nlNumber b(nlNumber x) {
  return (1/(L1*R*C));
}

nlNumber c(nlNumber x) {
  return ((1/L1)*((1/(L1*C))+(1/(L2*C))));
}

nlNumber d(nlNumber x) {
  return (1/(L1*L2*C)) * e(x);
}

nlNumber f2(nlNumber x, numberList* LastValues) {
  return d(x)-b(x)*nlGet(LastValues, 2)-c(x)*nlGet(LastValues, 1);
}

nlNumber printResult(nlNumber x, numberList* LastValues) {
  printf("%g,%g\n", x, nlGet(LastValues, 0));
}

main() {
  numberList* InitialValues;
  InitialValues = nlAllocList(3);
  nlSet(InitialValues, 0, 0);
  nlSet(InitialValues, 1, 0);
  nlSet(InitialValues, 2, 0);
  rungekutta(Start, End, StepSize, InitialValues, &f2, &printResult);
  nlFreeList(InitialValues);
}