from math import *

def iterative():
  I41 = 0.0
  d   = None
  while True:
    I41old = I41
    dold  = d

    t1  = sqrt(L*C) * asin(I41/k1)
    I23 = k2 * sin((ton-t1) / sqrt(L*C))
    t3  = sqrt(L*C) * asin(I23/k3)
    I41 = k4 * sin((toff-t3) / sqrt(L*C))

    d = I41 - I41old
    if dold != None and d*dold <= 0:
      break
    dold = d

  return t1, t3, I41, I23

def direct():
  r   = sqrt(L*C)
  s   = cos(ton / r)
  t   = cos(toff / r)
  u   = k1 / k4
  v   = k3 / k2
  uv  = u * v
  gu  = 1 - u**2
  gv  = 1 - v**2
  guv = 1 - uv**2
  p   = u*s - v*t
  f   = p*t*s + v*s - u*t
  b1  = 2 * v * (s - uv*t)
  b3  = 2 * u * (t - uv*s)
  a1  = guv * (gv + (v*t)**2 - s**2) + f*b1
  a3  = guv * (gu + (u*s)**2 - t**2) - f*b3
  q   = (1 - t**2) * (1 - s**2)
  d   = sqrt(q * (p**2 + gu*gv))
  e   = guv**2 - b1 * b3
  g1  = sqrt((a1 + b1*d) / e)
  g3  = sqrt((a3 + b3*d) / e)

  t1  = r * asin(g1)
  t3  = r * asin(g3)
  I41 = g1 * k1
  I23 = g3 * k3

  return t1, t3, I41, I23

# Beispiel

L    = 200e-6
C    = 100e-6
k1   = 2.5
k2   = 1
k3   = 3.2
k4   = 1
ton  = 70e-6
toff = 30e-6

for f in [iterative, direct]:
  print('t1 = %.9g  t3 = %.9g  I41 = %.9g  I23 = %.9g ' % f())
