diff --git a/docs/build/html/_modules/data.html b/docs/build/html/_modules/data.html index 679abdd..f93c332 100644 --- a/docs/build/html/_modules/data.html +++ b/docs/build/html/_modules/data.html @@ -89,7 +89,7 @@ object_data = pickle.load(input) # one load for every dump is needed to load all the data if verbose: print('found:') - print(object_data) + print(object_data) except IOError: object_data = None if verbose: diff --git a/docs/build/html/_modules/index.html b/docs/build/html/_modules/index.html index f359fb7..678a239 100644 --- a/docs/build/html/_modules/index.html +++ b/docs/build/html/_modules/index.html @@ -40,6 +40,7 @@
  • numerical.integration
  • numerical.ode
  • numerical.ode_model
    • +
  • time_of_day
    • diff --git a/docs/build/html/_modules/numerical/ode.html b/docs/build/html/_modules/numerical/ode.html index b827064..629bed9 100644 --- a/docs/build/html/_modules/numerical/ode.html +++ b/docs/build/html/_modules/numerical/ode.html @@ -37,8 +37,8 @@ # -*- coding: utf-8 -*- """Numerical solver of ordinary differential equations. -Solves the initial value problem for systems of first order ordinary differential -equations. +Solves the initial value problem for systems of first order +ordinary differential equations. :Date: 2015-09-21 @@ -68,7 +68,8 @@ :type x0: list :param t: time :type t: list - :param `*p`: parameters of the function (thickness, diameter, ...) + :param `*p`: parameters of the function (thickness, diameter, + ...) :param verbose: print information (default = False) :type verbose: bool @@ -83,14 +84,14 @@ .. math :: t_i = t_0 + i h ~,\quad i=1,2,\ldots,n - The derivative of the solution is approximated as the forward difference - equation + The derivative of the solution is approximated as the forward + difference equation .. math :: \dot{x}_i = f(t_i, x_i) = \frac{x_{i+1} - x_i}{t_{i+1}-t_i} - Therefore one step :math:`h` of the Euler method from :math:`t_i` to - :math:`t_{i+1}` is + Therefore one step :math:`h` of the Euler method from + :math:`t_i` to :math:`t_{i+1}` is .. math :: x_{i+1} &= x_i + (t_{i+1}-t_i) f(t_i, x_i) \\ @@ -119,7 +120,8 @@ .. math :: \dot{x} &= f(t,x) \\ \begin{bmatrix} \dot{x}_1 \\ \dot{x}_2 \end{bmatrix} &= - \begin{bmatrix} x_2 \\ m^{-1}(f(t) - d x_2 - k x_1) \end{bmatrix} \\ + \begin{bmatrix} x_2 \\ m^{-1}(f(t) - d x_2 - k x_1) + \end{bmatrix} \\ &= \begin{bmatrix} 0 \\ m^{-1} f(t) \end{bmatrix} + \begin{bmatrix} 0 & 1 \\ -m^{-1} k & -m^{-1} d \end{bmatrix} @@ -141,32 +143,39 @@ .. math :: \dot{x}_1 &= x_2 \\ - \dot{x}_2 &= m^{-1}(x_1)(f(t) - d(x_1,x_2) x_2 - k(x_1) x_1) \\ + \dot{x}_2 &= + m^{-1}(x_1)(f(t) - d(x_1,x_2) x_2 - k(x_1) x_1) \\ or .. math :: \dot{x} &= f(t,x) \\ \begin{bmatrix} \dot{x}_1 \\ \dot{x}_2 \end{bmatrix} &= - \begin{bmatrix} x_2 \\ m^{-1}(x_1)(f(t) - d(x_1,x_2) x_2 - k(x_1) x_1) \end{bmatrix} \\ + \begin{bmatrix} + x_2 \\ m^{-1}(x_1)(f(t) - d(x_1,x_2) x_2 - k(x_1) x_1) + \end{bmatrix} \\ &= \begin{bmatrix} 0 \\ m^{-1}(x_1) f(t) \end{bmatrix} + - \begin{bmatrix} 0 & 1 \\ -m^{-1}(x_1) k(x_1) & -m^{-1} d(x_1,x_2) \end{bmatrix} + \begin{bmatrix} + 0 & 1 \\ -m^{-1}(x_1) k(x_1) & -m^{-1} d(x_1,x_2) + \end{bmatrix} \begin{bmatrix} x_1 \\ x_2 \end{bmatrix} - The Euler method is a first-order method, - which means that the local error (error per step) is proportional to the - square of the step size, and the global error (error at a given time) is + The Euler method is a first-order method, which means that the + local error (error per step) is proportional to the square of + the step size, and the global error (error at a given time) is proportional to the step size. """ - x = zeros((len(t), len(x0))) # Preallocate array - x[0,:] = x0 # Initial condition gives solution at first t - for i in range(len(t)-1): # Calculation loop + x = zeros((len(t), len(x0))) # Preallocate array + x[0,:] = x0 # Initial condition gives solution at first t + for i in range(len(t)-1): # Calculation loop Dt = t[i+1]-t[i] dxdt = array(f(x[i,:], t[i], *p)) - x[i+1,:] = x[i,:] + dxdt*Dt # Approximate solution at next value of x + # Approximate solution at next value of x + x[i+1,:] = x[i,:] + dxdt*Dt if verbose: - print('Numerical integration of ODE using explicit first-order method (Euler / Runge-Kutta) was successful.') + print('Numerical integration of ODE using explicit ' + + 'first-order method (Euler / Runge-Kutta) was successful.') return x
    [docs]def e2(f, x0, t, *p, verbose=False): @@ -178,19 +187,22 @@ :type x0: list :param t: time :type t: list - :param `*p`: parameters of the function (thickness, diameter, ...) + :param `*p`: parameters of the function (thickness, diameter, + ...) :param verbose: print information (default = False) :type verbose: bool """ - x = zeros((len(t), len(x0))) # Preallocate array - x[0,:] = x0 # Initial condition gives solution at first t - for i in range(len(t)-1): # Calculation loop + x = zeros((len(t), len(x0))) # Preallocate array + x[0,:] = x0 # Initial condition gives solution at first t + for i in range(len(t)-1): # Calculation loop Dt = t[i+1]-t[i] k_1 = array(f(x[i,:], t[i], *p)) k_2 = array(f(x[i,:]+0.5*Dt*k_1, t[i]+0.5*Dt, *p)) - x[i+1,:] = x[i,:] + k_2*Dt # Approximate solution at next value of x + # Approximate solution at next value of x + x[i+1,:] = x[i,:] + k_2*Dt if verbose: - print('Numerical integration of ODE using explicit 2th-order method (Runge-Kutta) was successful.') + print('Numerical integration of ODE using explicit ' + + '2th-order method (Runge-Kutta) was successful.') return x
    [docs]def e4(f, x0, t, *p, verbose=False): @@ -202,7 +214,8 @@ :type x0: list :param t: time :type t: list - :param `*p`: parameters of the function (thickness, diameter, ...) + :param `*p`: parameters of the function (thickness, diameter, + ...) :param verbose: print information (default = False) :type verbose: bool """ @@ -214,70 +227,64 @@ k_2 = array(f(x[i,:]+0.5*Dt*k_1, t[i]+0.5*Dt, *p)) k_3 = array(f(x[i,:]+0.5*Dt*k_2, t[i]+0.5*Dt, *p)) k_4 = array(f(x[i,:]+k_3*Dt, t[i]+Dt, *p)) - x[i+1,:] = x[i,:] + 1./6*(k_1+2*k_2+2*k_3+k_4)*Dt # Approximate solution at next value of x + # Approximate solution at next value of x + x[i+1,:] = x[i,:] + 1./6*(k_1+2*k_2+2*k_3+k_4)*Dt if verbose: - print('Numerical integration of ODE using explicit 4th-order method (Runge-Kutta) was successful.') + print('Numerical integration of ODE using explicit ' + + '4th-order method (Runge-Kutta) was successful.') return x
    -
    [docs]def dxdt_Dt(f, x, t, Dt, *p): - r""" - :param f: :math:`f = \dot{x}` - :type f: function - :param Dt: :math:`\Delta{t}` +
    [docs]def fpi(f, xi, ti, ti1, *p, max_iterations=1000, tol=1e-9, + verbose=False): + r"""Fixed-point iteration. - :returns: :math:`\Delta x = \dot{x} \Delta t` - """ - return array(dxdt(x, t, *p)) * Dt
    - -
    [docs]def fixed_point_iteration(f, xi, t, max_iterations=1000, tol=1e-9, verbose=False): - r""" - :param f: the function to iterate :math:`f = \Delta{x}(t)` + :param f: the function to iterate :math:`f = \dot{x}(x,t)` :type f: function :param xi: initial condition :math:`x_i` :type xi: list - :param t: time :math:`t` - :type t: float - :param `*p`: parameters of the function (thickness, diameter, ...) + :param ti: time :math:`t_i` + :type ti: float + :param ti1: time :math:`t_{i+1}` + :type ti1: float + :param `*p`: parameters of the function (thickness, diameter, + ...) :param max_iterations: maximum number of iterations :type max_iterations: int - :param tol: tolerance against residuum (default = 1e-9) + :param tol: tolerance against residuum :math:`\varepsilon` + (default = 1e-9) :type tol: float :param verbose: print information (default = False) :type verbose: bool - :returns: :math:`x_{i+1}` + :returns: :math:`x_{i}` .. math :: - x_{i+1} = x_i + \Delta x - - .. seealso:: - :meth:`dxdt_Dt` for :math:`\Delta x` + x_{i,j=0} = x_{i} + + .. math :: + x_{i,j+1} = x_i + \dot{x}(x_{i,j}, t_{i+1})\cdot(t_{i+1}-t_i) + + .. math :: + \text{residuum} = \frac{\lVert x_{i,j+1}-x_{i,j}\rVert} + {\lVert x_{i,j+1} \rVert} < \varepsilon + + .. math :: + x_{i} = x_{i,j=\text{end}} """ - xi = x0 - for j in range(max_iterations): # Fixed-point iteration - Dx = array(f(xi, t, *p)) - xi1 = x0 + Dx # Approximate solution at next value of x - residuum = norm(xi1-xi)/norm(xi1) - xi = xi1 + xij = xi + for j in range(max_iterations): + dxdt = array(f(xij, ti1, *p)) + # Approximate solution at next value of x + xij1 = xi + dxdt * (ti1-ti) + residuum = norm(xij1-xij)/norm(xij1) + xij = xij1 if residuum < tol: break iterations = j+1 # number beginning with 1 therefore + 1 - return xi, iterations
    + return xij, iterations
    -
    [docs]def i1n(f, x0, t, *p, max_iterations=1000, tol=1e-9, verbose=False): - iterations = zeros((len(t), 1)) - x = zeros((len(t), len(x0))) # Preallocate array - x[0,:] = x0 # Initial condition gives solution at first t - for i in range(len(t)-1): - Dt = t[i+1]-t[i] - xi = x[i,:] - Dx = dxdt_Dt(f, xi, t[i+1], Dt, *p) - x[i+1,:], iterations[i] = fixed_point_iteration(Dx, xi, t, max_iterations, tol, verbose) - if verbose: - print('Numerical integration of ODE using implicite first-order method (Euler) was successful.') - return x, iterations
    - -
    [docs]def i1(f, x0, t, *p, max_iterations=1000, tol=1e-9, verbose=False): +
    [docs]def i1(f, x0, t, *p, max_iterations=1000, tol=1e-9, + verbose=False): r"""Implicite first-order method / backward Euler method. :param f: the function to solve @@ -286,7 +293,8 @@ :type x0: list :param t: time :type t: list - :param `*p`: parameters of the function (thickness, diameter, ...) + :param `*p`: parameters of the function (thickness, diameter, + ...) :param max_iterations: maximum number of iterations :type max_iterations: int :param tol: tolerance against residuum (default = 1e-9) @@ -298,26 +306,24 @@ """ iterations = zeros((len(t), 1)) x = zeros((len(t), len(x0))) # Preallocate array - x[0,:] = x0 # Initial condition gives solution at first t + x[0,:] = x0 # Initial condition gives solution at first t + # x(i+1) = x(i) + f(x(i+1), t(i+1)), exact value of + # f(x(i+1), t(i+1)) is not available therefore using + # Newton-Raphson method for i in range(len(t)-1): Dt = t[i+1]-t[i] xi = x[i,:] - # x(i+1) = x(i) + f(x(i+1), t(i+1)), exact value of f(x(i+1), t(i+1)) is not - # available therefor using Newton-Raphson method - for j in range(max_iterations): # Fixed-point iteration - dxdt = array(f(xi, t[i+1], *p)) - xi1 = x[i,:] + dxdt*Dt # Approximate solution at next value of x - residuum = norm(xi1-xi)/norm(xi1) - xi = xi1 - if residuum < tol: - break - iterations[i] = j+1 + xi, iteration = fpi(f, xi, t[i], t[i+1], *p, max_iterations, + tol, verbose) x[i+1,:] = xi + iterations[i] = iteration if verbose: - print('Numerical integration of ODE using implicite first-order method (Euler) was successful.') + print('Numerical integration of ODE using implicite ' + + 'first-order method (Euler) was successful.') return x, iterations
    -
    [docs]def newmark_newtonraphson(f, x0, xp0, xpp0, t, *p, gamma=.5, beta=.25, max_iterations=1000, tol=1e-9, verbose=False): +
    [docs]def newmark_newtonraphson(f, x0, xp0, xpp0, t, *p, gamma=.5, + beta=.25, max_iterations=1000, tol=1e-9, verbose=False): r"""Newmark method. :param f: the function to solve @@ -330,7 +336,8 @@ :type xpp0: list :param t: time :type t: list - :param `*p`: parameters of the function (thickness, diameter, ...) + :param `*p`: parameters of the function (thickness, diameter, + ...) :param gamma: newmark parameter for velocity (default = 0.5) :type gamma: float :param beta: newmark parameter for displacement (default = 0.25) @@ -346,9 +353,9 @@ x = zeros((len(t), len(x0))) # Preallocate array xp = zeros((len(t), len(xp0))) # Preallocate array xpp = zeros((len(t), len(xpp0))) # Preallocate array - x[0,:] = x0 # Initial condition gives solution at first t - xp[0,:] = xp0 # Initial condition gives solution at first t - xpp[0,:] = xpp0 # Initial condition gives solution at first t + x[0,:] = x0 # Initial condition gives solution at first t + xp[0,:] = xp0 # Initial condition gives solution at first t + xpp[0,:] = xpp0 # Initial condition gives solution at first t for i in range(len(t)-1): Dt = t[i+1]-t[i] @@ -361,7 +368,8 @@ j = 0 for j in range(max_iterations): # Fixed-point iteration #dxdt = array(f(t[i+1], x1, p)) - #x11 = x[i,:] + dxdt*Dt # Approximate solution at next value of x + # Approximate solution at next value of x + #x11 = x[i,:] + dxdt*Dt N, dN, dNp, dNpp = f(x1.reshape(-1,).tolist(), xp1.reshape(-1,).tolist(), xpp1.reshape(-1,).tolist(), @@ -370,7 +378,8 @@ print('divergiert') break - xpp11 = xpp1 - dot(inv(dNpp), (N + dot(dN, (x1-xi)) + dot(dNp, (xp1-xpi)))) + xpp11 = xpp1 - dot(inv(dNpp), (N + dot(dN, (x1-xi)) + \ + dot(dNp, (xp1-xpi)))) xp1 = xpi + Dt*( (1-gamma)*xppi + gamma*xpp11 ) x1 = xi + Dt*xpi + Dt**2*( (.5-beta)*xppi + beta*xpp11 ) @@ -384,11 +393,13 @@ xp[i+1,:] = xp1.reshape(-1,).tolist() x[i+1,:] = x1.reshape(-1,).tolist() if verbose: - print('Numerical integration of ODE using explicite newmark method was successful.') + print('Numerical integration of ODE using explicite ' + + 'newmark method was successful.') return x, xp, xpp, iterations
    # x = concatenate((x, xp, xpp), axis=1) -
    [docs]def newmark_newtonraphson_rdk(fnm, x0, xp0, xpp0, t, *p, gamma=.5, beta=.25, maxIterations=1000, tol=1e-9, verbose=False): +
    [docs]def newmark_newtonraphson_rdk(fnm, x0, xp0, xpp0, t, *p, gamma=.5, + beta=.25, max_iterations=1000, tol=1e-9, verbose=False): r"""Newmark method. :param f: the function to solve @@ -401,7 +412,8 @@ :type xpp0: list :param t: time :type t: list - :param `*p`: parameters of the function (thickness, diameter, ...) + :param `*p`: parameters of the function (thickness, diameter, + ...) :param gamma: newmark parameter for velocity (default = 0.5) :type gamma: float :param beta: newmark parameter for displacement (default = 0.25) @@ -417,13 +429,14 @@ x = zeros((len(t), len(x0))) # Preallocate array xp = zeros((len(t), len(xp0))) # Preallocate array xpp = zeros((len(t), len(xpp0))) # Preallocate array - x[0,:] = x0 # Initial condition gives solution at first t - xp[0,:] = xp0 # Initial condition gives solution at first t - xpp[0,:] = xpp0 # Initial condition gives solution at first t + x[0,:] = x0 # Initial condition gives solution at first t + xp[0,:] = xp0 # Initial condition gives solution at first t + xpp[0,:] = xpp0 # Initial condition gives solution at first t for i in range(len(t)-1): Dt = t[i+1]-t[i] - rm, rmx, rmxpp, rd, rdx, rdxp, rk, rkx, f = fnm(x[i,:], xp[i,:], xpp[i,:], t[i], *p) + rm, rmx, rmxpp, rd, rdx, rdxp, rk, rkx, f = fnm(x[i,:], + xp[i,:], xpp[i,:], t[i], *p) xi = x[i,:].reshape(3,1) xpi = xp[i,:].reshape(3,1) @@ -432,13 +445,16 @@ xp1 = xpi xpp1 = xppi j = 0 - for j in range(maxIterations): # Fixed-point iteration + for j in range(max_iterations): # Fixed-point iteration #dxdt = array(f(t[i+1], x1, p)) - #x11 = x[i,:] + dxdt*Dt # Approximate solution at next value of x + # Approximate solution at next value of x + #x11 = x[i,:] + dxdt*Dt r = (rmx+rdx+rkx)*Dt**2./4 + rdxp*Dt/2 + rmxpp - rp = f - (rm + dot(rmx, (Dt*xpi+Dt**2./4*xppi)) - dot(rmxpp, xppi) + \ - rd + dot(rdx, (Dt*xpi+Dt**2./4*xppi)) + dot(rdxp, Dt/2*xppi) + \ + rp = f - (rm + dot(rmx, (Dt*xpi+Dt**2./4*xppi)) - \ + dot(rmxpp, xppi) + \ + rd + dot(rdx, (Dt*xpi+Dt**2./4*xppi)) + \ + dot(rdxp, Dt/2*xppi) + \ rk + dot(rkx, (Dt*xpi+Dt**2./4*xppi)) ) xpp11 = dot(inv(r), rp) xp1 = xpi + Dt*( (1-gamma)*xppi + gamma*xpp11 ) @@ -454,7 +470,8 @@ xp[i+1,:] = xp1.reshape(-1,).tolist() x[i+1,:] = x1.reshape(-1,).tolist() if verbose: - print('Numerical integration of ODE using explicite newmark method was successful.') + print('Numerical integration of ODE using explicite ' + + 'newmark method was successful.') return x, xp, xpp, iterations
    # x = concatenate((x, xp, xpp), axis=1)
    diff --git a/docs/build/html/_modules/time_of_day.html b/docs/build/html/_modules/time_of_day.html new file mode 100644 index 0000000..b5afb74 --- /dev/null +++ b/docs/build/html/_modules/time_of_day.html @@ -0,0 +1,255 @@ + + + + + + + time_of_day — pylib 2019.5.19 documentation + + + + + + + + + + + + + + + + + + + + +
    +
    +
    + + +
    + +

    Source code for time_of_day

    +#!/usr/bin/env python
+# -*- coding: utf-8 -*-
+"""Calculate time.
+
+:Date: 2019-06-01
+
+.. module:: time_of_day
+  :platform: *nix, Windows
+  :synopsis: Calculate time.
+   
+.. moduleauthor:: Daniel Weschke <daniel.weschke@directbox.de>
+"""
+from __future__ import division, print_function, unicode_literals
+from time import struct_time, mktime
+
+
+
    [docs]def in_seconds(time):
+  """If time is `time.struct_time` convert to float seconds.
+
+  :param time: the time in seconds
+  :type time: float or `time.struct_time`
+  
+  :returns: the time in seconds
+  :rtype: float
+  """
+  if isinstance(time, struct_time):
+    time = mktime(time)
+  return time
    
+
+
    [docs]def seconds(time):
+  """The seconds of the time.
+
+  :param time: the time in seconds
+  :type time: float or `time.struct_time`
+  
+  :returns: seconds, range [0, 60]
+  :rtype: float
+  """
+  return in_seconds(time)%60
    
+
+
    [docs]def seconds_norm(time):
+  """The seconds normalized to 60 seconds.
+  
+  :param time: the time in seconds
+  :type time: float or `time.struct_time`
+  
+  :returns: the normalized seconds, range [0, 1]
+  :rtype: float
+  """
+  return seconds(time)/60
    
+
+
    [docs]def minutes(time):
+  """The minutes of the time.
+
+  :param time: the time in seconds
+  :type time: float or `time.struct_time`
+  
+  :returns: minutes, range [0, 60]
+  :rtype: float
+  """
+  return in_seconds(time)/60%60
    
+
+
    [docs]def minutes_norm(time):
+  """The minutes normalized to 60 minutes.
+  
+  :param time: the time in seconds
+  :type time: float or `time.struct_time`
+  
+  :returns: the normalized minutes, range [0, 1]
+  :rtype: float
+  """
+  return minutes(time)/60
    
+
+
    [docs]def hours(time):
+  """The hours of the time.
+
+  :param time: the time in seconds
+  :type time: float or `time.struct_time`
+  
+  :returns: hours, range [0, 24]
+  :rtype: float
+  """
+  return in_seconds(time)/60/60%24
    
+
+
    [docs]def hours_norm(time):
+  """The hours normalized to 24 hours.
+  
+  :param time: the time in seconds
+  :type time: float or `time.struct_time`
+  
+  :returns: the normalized hours, range [0, 1]
+  :rtype: float
+  """
+  return hours(time)/24
    
+
+
    [docs]def days(time):
+  """The days of the time (year).
+
+  :param time: the time in seconds
+  :type time: float or `time.struct_time`
+  
+  :returns: hours, range [0, 365.2425]
+  :rtype: float
+  """
+  return in_seconds(time)/60/60/24%365.2425
    
+
+
    [docs]def days_norm(time):
+  """The days normalized to 365.2425 (Gregorian, on average) days.
+  
+  :param time: the time in seconds
+  :type time: float or `time.struct_time`
+  
+  :returns: the normalized days, range [0, 1]
+  :rtype: float
+  """
+  return days(time)/365.2425
    
+
+
    [docs]def transform(time_norm, length, offset=0):
+  """Transform normalized time value to new length.
+  
+  :param position_norm: the normalized time value to transform
+  :type position_norm: float
+  :param length: the transformation
+  :type length: float
+  :param offset: the offset (default = 0)
+  :type offset: float
+  
+  :returns: the transformation value
+  :rtype: float
+  """
+  return time_norm*length + offset
    
+
+if __name__ == "__main__":
+  from time import time, gmtime, localtime
+  # time in seconds
+  t = time()
+  min = minutes(t)
+  h = hours(t)
+  min_norm = minutes_norm(t)
+  h_norm = hours_norm(t)
+  
+  print('min      ', min)
+  print('h        ', h)
+  print('min_norm ', min_norm)
+  print('h_norm   ', h_norm)
+  
+  x_len = 30
+  x_offset = -8
+  x_pos = transform(min_norm, x_len, x_offset)
+  print('m[-8,22] ', x_pos)
+  
+  y_len = 20
+  y_offset = -10
+  y_pos = transform(h_norm, y_len, y_offset)
+  print('h[-10,10]', y_pos)
+  
+
    + +
    + +
    +
    + +
    +
    +
    + ©2019, Daniel Weschke. + + | + Powered by Sphinx 2.0.1 + & Alabaster 0.7.12 + +
    + + + + + + \ No newline at end of file diff --git a/docs/build/html/_sources/modules.rst.txt b/docs/build/html/_sources/modules.rst.txt index 4679bea..f6d9421 100644 --- a/docs/build/html/_sources/modules.rst.txt +++ b/docs/build/html/_sources/modules.rst.txt @@ -8,3 +8,4 @@ src date geometry numerical + time_of_day diff --git a/docs/build/html/_sources/time_of_day.rst.txt b/docs/build/html/_sources/time_of_day.rst.txt new file mode 100644 index 0000000..4890e6b --- /dev/null +++ b/docs/build/html/_sources/time_of_day.rst.txt @@ -0,0 +1,7 @@ +time\_of\_day module +==================== + +.. automodule:: time_of_day + :members: + :undoc-members: + :show-inheritance: diff --git a/docs/build/html/genindex.html b/docs/build/html/genindex.html index 38c3bb8..bce3f3e 100644 --- a/docs/build/html/genindex.html +++ b/docs/build/html/genindex.html @@ -43,6 +43,7 @@ | E | F | G + | H | I | L | M @@ -84,18 +85,20 @@
  • data_read() (in module data)
  • data_store() (in module data) +
    • +
  • date (module), [1]
    • @@ -127,7 +130,7 @@ @@ -148,6 +151,18 @@ +

    H

    + + + +
    +

    I

      @@ -155,7 +170,7 @@
      -
    • i1n() (in module numerical.ode) +
    • in_seconds() (in module time_of_day)
    • integration (module)
      • @@ -174,6 +189,12 @@ +
      @@ -241,6 +262,12 @@

      S

      + @@ -249,10 +276,14 @@

      T

      diff --git a/docs/build/html/modules.html b/docs/build/html/modules.html index eaf6ac7..72af713 100644 --- a/docs/build/html/modules.html +++ b/docs/build/html/modules.html @@ -48,6 +48,7 @@
    • Module contents
      • +
    • time_of_day module
      • diff --git a/docs/build/html/numerical.html b/docs/build/html/numerical.html index 8e5ce87..e3b6b2b 100644 --- a/docs/build/html/numerical.html +++ b/docs/build/html/numerical.html @@ -197,30 +197,14 @@ b &= 1\end{split}\]

      numerical.ode module

      Numerical solver of ordinary differential equations.

      -

      Solves the initial value problem for systems of first order ordinary differential -equations.

      +

      Solves the initial value problem for systems of first order +ordinary differential equations.

      Date

      2015-09-21

      -
      -dxdt_Dt(f, x, t, Dt, *p)[source]
      -
      -
      Parameters
      -
        -

      • f (function) – \(f = \dot{x}\)

        • -

      • Dt – \(\Delta{t}\)

        • -
      -
      -
      Returns
      -

      \(\Delta x = \dot{x} \Delta t\)

      -
      -
      -
      - -
      e1(f, x0, t, *p, verbose=False)[source]

      Explicit first-order method / @@ -235,7 +219,8 @@ Euler-Cauchy-Verfahren / Euler-vorwärts-Verfahren

    • f (function) – the function to solve

    • x0 (list) – initial condition

    • t (list) – time

      • -

    • *p – parameters of the function (thickness, diameter, …)

      • +

    • *p – parameters of the function (thickness, diameter, +…)

    • verbose (bool) – print information (default = False)

      • @@ -247,12 +232,12 @@ x(t_0) &= x_0\end{split}\]

      Choose a value h for the size of every step and set

      \[t_i = t_0 + i h ~,\quad i=1,2,\ldots,n\]
      -

      The derivative of the solution is approximated as the forward difference -equation

      +

      The derivative of the solution is approximated as the forward +difference equation

      \[\dot{x}_i = f(t_i, x_i) = \frac{x_{i+1} - x_i}{t_{i+1}-t_i}\]
      -

      Therefore one step \(h\) of the Euler method from \(t_i\) to -\(t_{i+1}\) is

      +

      Therefore one step \(h\) of the Euler method from +\(t_i\) to \(t_{i+1}\) is

      \[\begin{split}x_{i+1} &= x_i + (t_{i+1}-t_i) f(t_i, x_i) \\ x_{i+1} &= x_i + h f(t_i, x_i) \\\end{split}\]
      @@ -272,7 +257,8 @@ x_2 &= \dot{u} &\quad \dot{x}_2 = \ddot{u} \\\end{split}\]
      \[\begin{split}\dot{x} &= f(t,x) \\ \begin{bmatrix} \dot{x}_1 \\ \dot{x}_2 \end{bmatrix} &= -\begin{bmatrix} x_2 \\ m^{-1}(f(t) - d x_2 - k x_1) \end{bmatrix} \\ +\begin{bmatrix} x_2 \\ m^{-1}(f(t) - d x_2 - k x_1) +\end{bmatrix} \\ &= \begin{bmatrix} 0 \\ m^{-1} f(t) \end{bmatrix} + \begin{bmatrix} 0 & 1 \\ -m^{-1} k & -m^{-1} d \end{bmatrix} @@ -288,19 +274,24 @@ x_2 &= \dot{u} &\quad \dot{x}_2 = \ddot{u} \\\end{split}\]

      becomes

      \[\begin{split}\dot{x}_1 &= x_2 \\ -\dot{x}_2 &= m^{-1}(x_1)(f(t) - d(x_1,x_2) x_2 - k(x_1) x_1) \\\end{split}\]
      +\dot{x}_2 &= + m^{-1}(x_1)(f(t) - d(x_1,x_2) x_2 - k(x_1) x_1) \\\end{split}\]

      or

      \[\begin{split}\dot{x} &= f(t,x) \\ \begin{bmatrix} \dot{x}_1 \\ \dot{x}_2 \end{bmatrix} &= -\begin{bmatrix} x_2 \\ m^{-1}(x_1)(f(t) - d(x_1,x_2) x_2 - k(x_1) x_1) \end{bmatrix} \\ +\begin{bmatrix} + x_2 \\ m^{-1}(x_1)(f(t) - d(x_1,x_2) x_2 - k(x_1) x_1) + \end{bmatrix} \\ &= \begin{bmatrix} 0 \\ m^{-1}(x_1) f(t) \end{bmatrix} + -\begin{bmatrix} 0 & 1 \\ -m^{-1}(x_1) k(x_1) & -m^{-1} d(x_1,x_2) \end{bmatrix} +\begin{bmatrix} + 0 & 1 \\ -m^{-1}(x_1) k(x_1) & -m^{-1} d(x_1,x_2) + \end{bmatrix} \begin{bmatrix} x_1 \\ x_2 \end{bmatrix}\end{split}\]
      -

      The Euler method is a first-order method, -which means that the local error (error per step) is proportional to the -square of the step size, and the global error (error at a given time) is +

      The Euler method is a first-order method, which means that the +local error (error per step) is proportional to the square of +the step size, and the global error (error at a given time) is proportional to the step size.

      @@ -314,7 +305,8 @@ proportional to the step size.

    • f (function) – the function to solve

    • x0 (list) – initial condition

    • t (list) – time

      • -

    • *p – parameters of the function (thickness, diameter, …)

      • +

    • *p – parameters of the function (thickness, diameter, +…)

    • verbose (bool) – print information (default = False)

      • @@ -331,7 +323,8 @@ proportional to the step size.

    • f (function) – the function to solve

    • x0 (list) – initial condition

    • t (list) – time

      • -

    • *p – parameters of the function (thickness, diameter, …)

      • +

    • *p – parameters of the function (thickness, diameter, +…)

    • verbose (bool) – print information (default = False)

      • @@ -339,30 +332,37 @@ proportional to the step size.

      -
      -fixed_point_iteration(f, xi, t, max_iterations=1000, tol=1e-09, verbose=False)[source]
      -
      +
      +fpi(f, xi, ti, ti1, *p, max_iterations=1000, tol=1e-09, verbose=False)[source]
      +

      Fixed-point iteration.

      +
      Parameters
        -

      • f (function) – the function to iterate \(f = \Delta{x}(t)\)

        • +

      • f (function) – the function to iterate \(f = \dot{x}(x,t)\)

      • xi (list) – initial condition \(x_i\)

        • -

      • t (float) – time \(t\)

        • -

      • *p – parameters of the function (thickness, diameter, …)

        • +

      • ti (float) – time \(t_i\)

        • +

      • ti1 (float) – time \(t_{i+1}\)

        • +

      • *p – parameters of the function (thickness, diameter, +…)

      • max_iterations (int) – maximum number of iterations

        • -

      • tol (float) – tolerance against residuum (default = 1e-9)

        • +

      • tol (float) – tolerance against residuum \(\varepsilon\) +(default = 1e-9)

      • verbose (bool) – print information (default = False)

      Returns
      -

      \(x_{i+1}\)

      +

      \(x_{i}\)

      -\[x_{i+1} = x_i + \Delta x\]
      -
      -

      See also

      -

      dxdt_Dt() for \(\Delta x\)

      -
      +\[x_{i,j=0} = x_{i}\] +
      +\[x_{i,j+1} = x_i + \dot{x}(x_{i,j}, t_{i+1})\cdot(t_{i+1}-t_i)\]
      +
      +\[\text{residuum} = \frac{\lVert x_{i,j+1}-x_{i,j}\rVert} + {\lVert x_{i,j+1} \rVert} < \varepsilon\]
      +
      +\[x_{i} = x_{i,j=\text{end}}\]
      @@ -375,7 +375,8 @@ proportional to the step size.

    • f (function) – the function to solve

    • x0 (list) – initial condition

    • t (list) – time

      • -

    • *p – parameters of the function (thickness, diameter, …)

      • +

    • *p – parameters of the function (thickness, diameter, +…)

    • max_iterations (int) – maximum number of iterations

    • tol (float) – tolerance against residuum (default = 1e-9)

    • verbose (bool) – print information (default = False)

      • @@ -385,11 +386,6 @@ proportional to the step size.

      The backward Euler method has order one and is A-stable.

      -
      -
      -i1n(f, x0, t, *p, max_iterations=1000, tol=1e-09, verbose=False)[source]
      -
      -
      newmark_newtonraphson(f, x0, xp0, xpp0, t, *p, gamma=0.5, beta=0.25, max_iterations=1000, tol=1e-09, verbose=False)[source]
      @@ -402,7 +398,8 @@ proportional to the step size.

    • xp0 (list) – initial condition

    • xpp0 (list) – initial condition

    • t (list) – time

      • -

    • *p – parameters of the function (thickness, diameter, …)

      • +

    • *p – parameters of the function (thickness, diameter, +…)

    • gamma (float) – newmark parameter for velocity (default = 0.5)

    • beta (float) – newmark parameter for displacement (default = 0.25)

    • max_iterations (int) – maximum number of iterations

      • @@ -415,7 +412,7 @@ proportional to the step size.

      -newmark_newtonraphson_rdk(fnm, x0, xp0, xpp0, t, *p, gamma=0.5, beta=0.25, maxIterations=1000, tol=1e-09, verbose=False)[source]
      +newmark_newtonraphson_rdk(fnm, x0, xp0, xpp0, t, *p, gamma=0.5, beta=0.25, max_iterations=1000, tol=1e-09, verbose=False)[source]

      Newmark method.

      Parameters
      @@ -425,7 +422,8 @@ proportional to the step size.

    • xp0 (list) – initial condition

    • xpp0 (list) – initial condition

    • t (list) – time

      • -

    • *p – parameters of the function (thickness, diameter, …)

      • +

    • *p – parameters of the function (thickness, diameter, +…)

    • gamma (float) – newmark parameter for velocity (default = 0.5)

    • beta (float) – newmark parameter for displacement (default = 0.25)

    • max_iterations (int) – maximum number of iterations

      • diff --git a/docs/build/html/objects.inv b/docs/build/html/objects.inv index 710740c..b6e0d3d 100644 --- a/docs/build/html/objects.inv +++ b/docs/build/html/objects.inv @@ -2,6 +2,6 @@ # Project: pylib # Version: # The remainder of this file is compressed using zlib. -xڭ0<h! 40CAQcGD# A&< ddH4}fL;mŬ&"e2uHZ*Dp/!LQQmV3inQpY $dfIfv^f{mD%"k:o~7 {i%Rj![꼱p[[*X*8 Np K0Jl6X"w}aߎ@jPȳZE`%gꢘ   AU +KRN4H._o{*aglklA o:V(4My;vx'>?AM 3V/G,.Sy n @t1e Ҩ)IOTsM G]mxFm)HbW35* \ No newline at end of file diff --git a/docs/build/html/py-modindex.html b/docs/build/html/py-modindex.html index df9bb3e..3a2d315 100644 --- a/docs/build/html/py-modindex.html +++ b/docs/build/html/py-modindex.html @@ -48,7 +48,8 @@ g | i | n | - o + o | + t
      @@ -137,6 +138,14 @@ + + + + +
      ode_model (*nix, Windows) Models of ordinary differential equations.
       
      + t
      + time_of_day (*nix, Windows) + Calculate time.
      diff --git a/docs/build/html/searchindex.js b/docs/build/html/searchindex.js index c9c5aa6..fa5c2a4 100644 --- a/docs/build/html/searchindex.js +++ b/docs/build/html/searchindex.js @@ -1 +1 @@ -Search.setIndex({docnames:["data","date","geometry","index","modules","numerical"],envversion:{"sphinx.domains.c":1,"sphinx.domains.changeset":1,"sphinx.domains.cpp":1,"sphinx.domains.javascript":1,"sphinx.domains.math":2,"sphinx.domains.python":1,"sphinx.domains.rst":1,"sphinx.domains.std":1,"sphinx.ext.viewcode":1,sphinx:56},filenames:["data.rst","date.rst","geometry.rst","index.rst","modules.rst","numerical.rst"],objects:{"":{data:[0,0,0,"-"],date:[1,0,0,"-"],fit:[5,0,0,"-"],geometry:[2,0,0,"-"],integration:[5,0,0,"-"],numerical:[5,0,0,"-"],ode:[5,0,0,"-"],ode_model:[5,0,0,"-"]},"numerical.fit":{gauss:[5,1,1,""],gauss_fit:[5,1,1,""]},"numerical.integration":{trapez:[5,1,1,""]},"numerical.ode":{dxdt_Dt:[5,1,1,""],e1:[5,1,1,""],e2:[5,1,1,""],e4:[5,1,1,""],fixed_point_iteration:[5,1,1,""],i1:[5,1,1,""],i1n:[5,1,1,""],newmark_newtonraphson:[5,1,1,""],newmark_newtonraphson_rdk:[5,1,1,""]},"numerical.ode_model":{disk:[5,1,1,""],disk_nm:[5,1,1,""],disk_nmmdk:[5,1,1,""]},data:{data_load:[0,1,1,""],data_read:[0,1,1,""],data_store:[0,1,1,""],main:[0,1,1,""]},date:{"gau\u00dfsche_osterformel":[1,1,1,""],ascension_of_jesus:[1,1,1,""],easter_friday:[1,1,1,""],easter_monday:[1,1,1,""],easter_sunday:[1,1,1,""],pentecost:[1,1,1,""]},geometry:{cubic:[2,1,1,""],cubic_deg:[2,1,1,""],line:[2,1,1,""],points:[2,1,1,""],rectangle:[2,1,1,""],rotate:[2,1,1,""],rotate_deg:[2,1,1,""],square:[2,1,1,""],translate:[2,1,1,""]},numerical:{fit:[5,0,0,"-"],integration:[5,0,0,"-"],ode:[5,0,0,"-"],ode_model:[5,0,0,"-"]}},objnames:{"0":["py","module","Python module"],"1":["py","function","Python function"]},objtypes:{"0":"py:module","1":"py:function"},terms:{"1st":5,"2nd":5,"4th":5,"9fsche_osterformel":1,"case":5,"default":[0,5],"f\u00fcr":1,"float":[2,5],"fr\u00fchling":1,"function":[0,5],"gau\u00dfsch":1,"gau\u00dfsche_osterformel":1,"int":[0,1,2,5],"korrekturgr\u00f6\u00df":1,"m\u00e4rz":1,"m\u00e4rzdatum":1,"return":[0,1,2,5],"s\u00e4kular":1,"s\u00e4kularzahl":1,"vorw\u00e4rt":5,Das:1,The:[2,5],against:5,algorithmu:1,als:1,amplitud:5,analyt:5,angl:2,angle1:2,angle2:2,approx:5,approxim:5,april:1,area:5,around:2,ascens:1,ascension_of_jesu:1,ascii:0,backward:5,base:5,becom:5,begin:5,beta:5,binari:0,bmatrix:5,bool:[0,5],calcul:[1,5],cauchi:5,center:2,choos:5,column:0,composit:5,condit:5,content:4,counterclockwis:2,cubic:2,cubic_deg:2,dai:1,data:[4,5],data_load:0,data_read:0,data_stor:0,date:[2,4,5],datetim:1,datum:1,ddot:5,definit:5,degre:2,delta:5,den:1,der:1,deriv:5,des:1,describ:5,deviat:5,diamet:5,die:1,differ:5,differenti:5,disk:5,disk_nm:5,disk_nmmdk:5,displac:5,distribut:5,divid:5,dot:5,dxdt_dt:5,easter:1,easter_fridai:1,easter_mondai:1,easter_sundai:1,eccentr:5,edg:2,end:[2,5],entfernung:1,equal:5,equalii:5,equat:5,error:5,ersten:1,euler:5,everi:5,exampl:5,expect:5,explicit:5,explizit:5,f_n:5,fals:[0,5],fassregel:5,file:0,file_nam:0,filenam:0,first:[0,5],fit:4,fixed_point_iter:5,float64:5,fnm:5,foral:5,forward:5,fourth:5,frac:[2,5],fridai:1,from:[0,5],full:5,fwhm:5,gamma:5,gau:1,gauss:5,gauss_fit:5,geometri:4,gilt:1,given:[2,5],global:5,govern:5,gregorianischen:1,half:5,has:5,height:2,http:1,i1n:5,implicit:5,index:[0,3],inform:[0,5],initi:5,integr:4,integrand:5,interv:5,iter:5,jahr:1,jesu:1,johann:5,kalend:1,kalendarisch:1,keim:1,kepler:5,keplersch:5,kutta:5,kwarg:2,lambda:5,ldot:5,left:5,leftrightarrow:5,leq:5,limit:5,limits_:5,limits_a:5,line:2,list:[0,5],load:0,local:5,lower:5,main:0,march:1,mathmat:5,mathrm:5,max_iter:5,maximum:5,maxiter:5,mean:5,method:5,model:5,modul:[3,4],mondai:1,mondparamet:1,mondschaltung:1,ndarrai:5,necessari:5,necessarili:5,newmark:5,newmark_newtonraphson:5,newmark_newtonraphson_rdk:5,none:5,number:[2,5],numer:4,numerisch:5,numpi:5,object:[0,2],object_data:0,ode:4,ode_model:4,offset:5,one:[2,5],option:2,order:5,ordinari:5,org:1,origin:2,osterentfernung:1,osterformel:1,ostergrenz:1,ostersonntag:1,other:2,packag:4,page:3,paramet:[0,1,2,5],pentecost:1,per:5,point1:2,point1_i:2,point1_x:2,point2:2,point2_i:2,point2_x:2,point3:2,point4:2,point:[2,5],point_i:2,point_x:2,points1_i:2,polygon:2,polygonzugverfahren:5,popt:5,posit:5,print:5,problem:5,program:0,proport:5,pts:2,quad:5,quadratur:5,radian:2,read:0,rearrang:2,rectangl:2,residuum:5,right:5,rotat:[2,5],rotate_deg:2,rule:5,rung:5,sampl:2,sample_point1_x:2,sample_point2_i:2,sample_point2_x:2,sample_points1_i:2,save_valu:5,sche:5,search:3,second:5,segment:5,set:5,should:2,side:5,sigma:5,simpson:5,simpsonregel:5,size:[2,5],slope:2,solut:5,solv:5,solver:5,sonnenschaltung:1,sonntag:1,sourc:[0,1,2,5],space:5,spacial:1,special:5,specifi:5,sqrt:5,squar:[2,5],stabl:5,standard:5,step:5,store:0,str:0,subinterv:5,submodul:4,sum:5,sundai:1,system:5,t_0:5,t_i:5,tagen:1,therefor:5,thick:5,thoma:5,time:5,tol:5,toler:5,torqu:5,translat:2,trapez:5,trapezium:5,trapezoid:5,trapezregel:5,tupl:[0,2,5],two:5,type:[0,1,2,5],upper:5,used:5,using:5,usw:1,valu:5,variabl:[1,5],varianc:5,vec:2,vector:2,veloc:5,verbos:[0,5],verfahren:5,vert_0:5,vert_a:5,vertic:5,vollmond:1,von:1,where:5,which:5,width:[2,5],wiki:1,wikipedia:1,write:0,x_0:5,x_1:[2,5],x_2:[2,5],x_column:0,x_fit:5,x_i:5,x_n:5,xp0:5,xpn:5,xpp0:5,xppn:5,y_1:2,y_2:2,y_column:0,y_fit:5,year:1},titles:["data module","date module","geometry module","Welcome to pylib\u2019s documentation!","src","numerical package"],titleterms:{content:5,data:0,date:1,document:3,fit:5,geometri:2,indic:3,integr:5,modul:[0,1,2,5],numer:5,ode:5,ode_model:5,packag:5,pylib:3,src:4,submodul:5,tabl:3,welcom:3}}) \ No newline at end of file +Search.setIndex({docnames:["data","date","geometry","index","modules","numerical","time_of_day"],envversion:{"sphinx.domains.c":1,"sphinx.domains.changeset":1,"sphinx.domains.cpp":1,"sphinx.domains.javascript":1,"sphinx.domains.math":2,"sphinx.domains.python":1,"sphinx.domains.rst":1,"sphinx.domains.std":1,"sphinx.ext.viewcode":1,sphinx:56},filenames:["data.rst","date.rst","geometry.rst","index.rst","modules.rst","numerical.rst","time_of_day.rst"],objects:{"":{data:[0,0,0,"-"],date:[1,0,0,"-"],fit:[5,0,0,"-"],geometry:[2,0,0,"-"],integration:[5,0,0,"-"],numerical:[5,0,0,"-"],ode:[5,0,0,"-"],ode_model:[5,0,0,"-"],time_of_day:[6,0,0,"-"]},"numerical.fit":{gauss:[5,1,1,""],gauss_fit:[5,1,1,""]},"numerical.integration":{trapez:[5,1,1,""]},"numerical.ode":{e1:[5,1,1,""],e2:[5,1,1,""],e4:[5,1,1,""],fpi:[5,1,1,""],i1:[5,1,1,""],newmark_newtonraphson:[5,1,1,""],newmark_newtonraphson_rdk:[5,1,1,""]},"numerical.ode_model":{disk:[5,1,1,""],disk_nm:[5,1,1,""],disk_nmmdk:[5,1,1,""]},data:{data_load:[0,1,1,""],data_read:[0,1,1,""],data_store:[0,1,1,""],main:[0,1,1,""]},date:{"gau\u00dfsche_osterformel":[1,1,1,""],ascension_of_jesus:[1,1,1,""],easter_friday:[1,1,1,""],easter_monday:[1,1,1,""],easter_sunday:[1,1,1,""],pentecost:[1,1,1,""]},geometry:{cubic:[2,1,1,""],cubic_deg:[2,1,1,""],line:[2,1,1,""],points:[2,1,1,""],rectangle:[2,1,1,""],rotate:[2,1,1,""],rotate_deg:[2,1,1,""],square:[2,1,1,""],translate:[2,1,1,""]},numerical:{fit:[5,0,0,"-"],integration:[5,0,0,"-"],ode:[5,0,0,"-"],ode_model:[5,0,0,"-"]},time_of_day:{days:[6,1,1,""],days_norm:[6,1,1,""],hours:[6,1,1,""],hours_norm:[6,1,1,""],in_seconds:[6,1,1,""],minutes:[6,1,1,""],minutes_norm:[6,1,1,""],seconds:[6,1,1,""],seconds_norm:[6,1,1,""],transform:[6,1,1,""]}},objnames:{"0":["py","module","Python module"],"1":["py","function","Python function"]},objtypes:{"0":"py:module","1":"py:function"},terms:{"1st":5,"2nd":5,"4th":5,"9fsche_osterformel":1,"case":5,"default":[0,5,6],"f\u00fcr":1,"float":[2,5,6],"fr\u00fchling":1,"function":[0,5],"gau\u00dfsch":1,"gau\u00dfsche_osterformel":1,"int":[0,1,2,5],"korrekturgr\u00f6\u00df":1,"m\u00e4rz":1,"m\u00e4rzdatum":1,"new":6,"return":[0,1,2,5,6],"s\u00e4kular":1,"s\u00e4kularzahl":1,"vorw\u00e4rt":5,Das:1,The:[2,5,6],against:5,algorithmu:1,als:1,amplitud:5,analyt:5,angl:2,angle1:2,angle2:2,approx:5,approxim:5,april:1,area:5,around:2,ascens:1,ascension_of_jesu:1,ascii:0,averag:6,backward:5,base:5,becom:5,begin:5,beta:5,binari:0,bmatrix:5,bool:[0,5],calcul:[1,5,6],cauchi:5,cdot:5,center:2,choos:5,column:0,composit:5,condit:5,content:4,convert:6,counterclockwis:2,cubic:2,cubic_deg:2,dai:[1,6],data:[4,5],data_load:0,data_read:0,data_stor:0,date:[2,4,5,6],datetim:1,datum:1,days_norm:6,ddot:5,definit:5,degre:2,delta:[],den:1,der:1,deriv:5,des:1,describ:5,deviat:5,diamet:5,die:1,differ:5,differenti:5,disk:5,disk_nm:5,disk_nmmdk:5,displac:5,distribut:5,divid:5,dot:5,dxdt_dt:[],easter:1,easter_fridai:1,easter_mondai:1,easter_sundai:1,eccentr:5,edg:2,end:[2,5],entfernung:1,equal:5,equalii:5,equat:5,error:5,ersten:1,euler:5,everi:5,exampl:5,expect:5,explicit:5,explizit:5,f_n:5,fals:[0,5],fassregel:5,file:0,file_nam:0,filenam:0,first:[0,5],fit:4,fix:5,fixed_point_iter:[],float64:5,fnm:5,foral:5,forward:5,fourth:5,fpi:5,frac:[2,5],fridai:1,from:[0,5],full:5,fwhm:5,gamma:5,gau:1,gauss:5,gauss_fit:5,geometri:4,gilt:1,given:[2,5],global:5,govern:5,gregorian:6,gregorianischen:1,half:5,has:5,height:2,hour:6,hours_norm:6,http:1,i1n:[],i1new:[],implicit:5,in_second:6,index:[0,3],inform:[0,5],initi:5,integr:4,integrand:5,interv:5,iter:5,jahr:1,jesu:1,johann:5,kalend:1,kalendarisch:1,keim:1,kepler:5,keplersch:5,kutta:5,kwarg:2,lambda:5,ldot:5,left:5,leftrightarrow:5,length:6,leq:5,limit:5,limits_:5,limits_a:5,line:2,list:[0,5],load:0,local:5,lower:5,lvert:5,main:0,march:1,mathmat:5,mathrm:5,max_iter:5,maximum:5,maxiter:[],mean:5,method:5,minut:6,minutes_norm:6,model:5,modul:[3,4],mondai:1,mondparamet:1,mondschaltung:1,ndarrai:5,necessari:5,necessarili:5,newmark:5,newmark_newtonraphson:5,newmark_newtonraphson_rdk:5,none:5,norm:[],normal:6,number:[2,5],numer:4,numerisch:5,numpi:5,object:[0,2],object_data:0,ode:4,ode_model:4,offset:[5,6],one:[2,5],option:2,order:5,ordinari:5,org:1,origin:2,osterentfernung:1,osterformel:1,ostergrenz:1,ostersonntag:1,other:2,packag:4,page:3,paramet:[0,1,2,5,6],pentecost:1,per:5,point1:2,point1_i:2,point1_x:2,point2:2,point2_i:2,point2_x:2,point3:2,point4:2,point:[2,5],point_i:2,point_x:2,points1_i:2,polygon:2,polygonzugverfahren:5,popt:5,posit:5,position_norm:6,print:5,problem:5,program:0,proport:5,pts:2,quad:5,quadratur:5,radian:2,rang:6,read:0,rearrang:2,rectangl:2,residuum:5,right:5,rotat:[2,5],rotate_deg:2,rule:5,rung:5,rvert:5,sampl:2,sample_point1_x:2,sample_point2_i:2,sample_point2_x:2,sample_points1_i:2,save_valu:5,sche:5,search:3,second:[5,6],seconds_norm:6,segment:5,set:5,should:2,side:5,sigma:5,simpson:5,simpsonregel:5,size:[2,5],slope:2,solut:5,solv:5,solver:5,sonnenschaltung:1,sonntag:1,sourc:[0,1,2,5,6],space:5,spacial:1,special:5,specifi:5,sqrt:5,squar:[2,5],stabl:5,standard:5,step:5,store:0,str:0,struct_tim:6,subinterv:5,submodul:4,sum:5,sundai:1,system:5,t_0:5,t_i:5,tagen:1,text:5,therefor:5,thick:5,thoma:5,ti1:5,time:[5,6],time_norm:6,time_of_dai:4,tol:5,toler:5,torqu:5,transform:6,translat:2,trapez:5,trapezium:5,trapezoid:5,trapezregel:5,tupl:[0,2,5],two:5,type:[0,1,2,5,6],upper:5,used:5,using:5,usw:1,valu:[5,6],varepsilon:5,variabl:[1,5],varianc:5,vec:2,vector:2,veloc:5,verbos:[0,5],verfahren:5,vert:[],vert_0:5,vert_a:5,vertic:5,vollmond:1,von:1,where:5,which:5,width:[2,5],wiki:1,wikipedia:1,write:0,x_0:5,x_1:[2,5],x_2:[2,5],x_column:0,x_fit:5,x_i:5,x_n:5,xp0:5,xpn:5,xpp0:5,xppn:5,y_1:2,y_2:2,y_column:0,y_fit:5,year:[1,6]},titles:["data module","date module","geometry module","Welcome to pylib\u2019s documentation!","src","numerical package","time_of_day module"],titleterms:{content:5,data:0,date:1,document:3,fit:5,geometri:2,indic:3,integr:5,modul:[0,1,2,5,6],numer:5,ode:5,ode_model:5,packag:5,pylib:3,src:4,submodul:5,tabl:3,time_of_dai:6,welcom:3}}) \ No newline at end of file diff --git a/docs/build/html/time_of_day.html b/docs/build/html/time_of_day.html new file mode 100644 index 0000000..12c9f53 --- /dev/null +++ b/docs/build/html/time_of_day.html @@ -0,0 +1,282 @@ + + + + + + + time_of_day module — pylib 2019.5.19 documentation + + + + + + + + + + + + + + + + + + + + +
      +
      +
      + + +
      + +
      +

      time_of_day module

      +

      Calculate time.

      +
      +
      Date
      +

      2019-06-01

      +
      +
      +
      +
      +days(time)[source]
      +

      The days of the time (year).

      +
      +
      Parameters
      +

      time (float or time.struct_time) – the time in seconds

      +
      +
      Returns
      +

      hours, range [0, 365.2425]

      +
      +
      Return type
      +

      float

      +
      +
      +
      + +
      +
      +days_norm(time)[source]
      +

      The days normalized to 365.2425 (Gregorian, on average) days.

      +
      +
      Parameters
      +

      time (float or time.struct_time) – the time in seconds

      +
      +
      Returns
      +

      the normalized days, range [0, 1]

      +
      +
      Return type
      +

      float

      +
      +
      +
      + +
      +
      +hours(time)[source]
      +

      The hours of the time.

      +
      +
      Parameters
      +

      time (float or time.struct_time) – the time in seconds

      +
      +
      Returns
      +

      hours, range [0, 24]

      +
      +
      Return type
      +

      float

      +
      +
      +
      + +
      +
      +hours_norm(time)[source]
      +

      The hours normalized to 24 hours.

      +
      +
      Parameters
      +

      time (float or time.struct_time) – the time in seconds

      +
      +
      Returns
      +

      the normalized hours, range [0, 1]

      +
      +
      Return type
      +

      float

      +
      +
      +
      + +
      +
      +in_seconds(time)[source]
      +

      If time is time.struct_time convert to float seconds.

      +
      +
      Parameters
      +

      time (float or time.struct_time) – the time in seconds

      +
      +
      Returns
      +

      the time in seconds

      +
      +
      Return type
      +

      float

      +
      +
      +
      + +
      +
      +minutes(time)[source]
      +

      The minutes of the time.

      +
      +
      Parameters
      +

      time (float or time.struct_time) – the time in seconds

      +
      +
      Returns
      +

      minutes, range [0, 60]

      +
      +
      Return type
      +

      float

      +
      +
      +
      + +
      +
      +minutes_norm(time)[source]
      +

      The minutes normalized to 60 minutes.

      +
      +
      Parameters
      +

      time (float or time.struct_time) – the time in seconds

      +
      +
      Returns
      +

      the normalized minutes, range [0, 1]

      +
      +
      Return type
      +

      float

      +
      +
      +
      + +
      +
      +seconds(time)[source]
      +

      The seconds of the time.

      +
      +
      Parameters
      +

      time (float or time.struct_time) – the time in seconds

      +
      +
      Returns
      +

      seconds, range [0, 60]

      +
      +
      Return type
      +

      float

      +
      +
      +
      + +
      +
      +seconds_norm(time)[source]
      +

      The seconds normalized to 60 seconds.

      +
      +
      Parameters
      +

      time (float or time.struct_time) – the time in seconds

      +
      +
      Returns
      +

      the normalized seconds, range [0, 1]

      +
      +
      Return type
      +

      float

      +
      +
      +
      + +
      +
      +transform(time_norm, length, offset=0)[source]
      +

      Transform normalized time value to new length.

      +
      +
      Parameters
      +
        +

      • position_norm (float) – the normalized time value to transform

        • +

      • length (float) – the transformation

        • +

      • offset (float) – the offset (default = 0)

        • +
      +
      +
      Returns
      +

      the transformation value

      +
      +
      Return type
      +

      float

      +
      +
      +
      + +
      + + +
      + +
      +
      + +
      +
      +
      + ©2019, Daniel Weschke. + + | + Powered by Sphinx 2.0.1 + & Alabaster 0.7.12 + + | + Page source +
      + + + + + + \ No newline at end of file diff --git a/docs/source/modules.rst b/docs/source/modules.rst index 4679bea..f6d9421 100644 --- a/docs/source/modules.rst +++ b/docs/source/modules.rst @@ -8,3 +8,4 @@ src date geometry numerical + time_of_day diff --git a/docs/source/time_of_day.rst b/docs/source/time_of_day.rst new file mode 100644 index 0000000..4890e6b --- /dev/null +++ b/docs/source/time_of_day.rst @@ -0,0 +1,7 @@ +time\_of\_day module +==================== + +.. automodule:: time_of_day + :members: + :undoc-members: + :show-inheritance: diff --git a/src/numerical/ode.py b/src/numerical/ode.py index 44a0caf..2673f9b 100644 --- a/src/numerical/ode.py +++ b/src/numerical/ode.py @@ -2,8 +2,8 @@ # -*- coding: utf-8 -*- """Numerical solver of ordinary differential equations. -Solves the initial value problem for systems of first order ordinary differential -equations. +Solves the initial value problem for systems of first order +ordinary differential equations. :Date: 2015-09-21 @@ -33,7 +33,8 @@ def e1(f, x0, t, *p, verbose=False): :type x0: list :param t: time :type t: list - :param `*p`: parameters of the function (thickness, diameter, ...) + :param `*p`: parameters of the function (thickness, diameter, + ...) :param verbose: print information (default = False) :type verbose: bool @@ -48,14 +49,14 @@ def e1(f, x0, t, *p, verbose=False): .. math :: t_i = t_0 + i h ~,\quad i=1,2,\ldots,n - The derivative of the solution is approximated as the forward difference - equation + The derivative of the solution is approximated as the forward + difference equation .. math :: \dot{x}_i = f(t_i, x_i) = \frac{x_{i+1} - x_i}{t_{i+1}-t_i} - Therefore one step :math:`h` of the Euler method from :math:`t_i` to - :math:`t_{i+1}` is + Therefore one step :math:`h` of the Euler method from + :math:`t_i` to :math:`t_{i+1}` is .. math :: x_{i+1} &= x_i + (t_{i+1}-t_i) f(t_i, x_i) \\ @@ -84,7 +85,8 @@ def e1(f, x0, t, *p, verbose=False): .. math :: \dot{x} &= f(t,x) \\ \begin{bmatrix} \dot{x}_1 \\ \dot{x}_2 \end{bmatrix} &= - \begin{bmatrix} x_2 \\ m^{-1}(f(t) - d x_2 - k x_1) \end{bmatrix} \\ + \begin{bmatrix} x_2 \\ m^{-1}(f(t) - d x_2 - k x_1) + \end{bmatrix} \\ &= \begin{bmatrix} 0 \\ m^{-1} f(t) \end{bmatrix} + \begin{bmatrix} 0 & 1 \\ -m^{-1} k & -m^{-1} d \end{bmatrix} @@ -106,32 +108,39 @@ def e1(f, x0, t, *p, verbose=False): .. math :: \dot{x}_1 &= x_2 \\ - \dot{x}_2 &= m^{-1}(x_1)(f(t) - d(x_1,x_2) x_2 - k(x_1) x_1) \\ + \dot{x}_2 &= + m^{-1}(x_1)(f(t) - d(x_1,x_2) x_2 - k(x_1) x_1) \\ or .. math :: \dot{x} &= f(t,x) \\ \begin{bmatrix} \dot{x}_1 \\ \dot{x}_2 \end{bmatrix} &= - \begin{bmatrix} x_2 \\ m^{-1}(x_1)(f(t) - d(x_1,x_2) x_2 - k(x_1) x_1) \end{bmatrix} \\ + \begin{bmatrix} + x_2 \\ m^{-1}(x_1)(f(t) - d(x_1,x_2) x_2 - k(x_1) x_1) + \end{bmatrix} \\ &= \begin{bmatrix} 0 \\ m^{-1}(x_1) f(t) \end{bmatrix} + - \begin{bmatrix} 0 & 1 \\ -m^{-1}(x_1) k(x_1) & -m^{-1} d(x_1,x_2) \end{bmatrix} + \begin{bmatrix} + 0 & 1 \\ -m^{-1}(x_1) k(x_1) & -m^{-1} d(x_1,x_2) + \end{bmatrix} \begin{bmatrix} x_1 \\ x_2 \end{bmatrix} - The Euler method is a first-order method, - which means that the local error (error per step) is proportional to the - square of the step size, and the global error (error at a given time) is + The Euler method is a first-order method, which means that the + local error (error per step) is proportional to the square of + the step size, and the global error (error at a given time) is proportional to the step size. """ - x = zeros((len(t), len(x0))) # Preallocate array - x[0,:] = x0 # Initial condition gives solution at first t - for i in range(len(t)-1): # Calculation loop + x = zeros((len(t), len(x0))) # Preallocate array + x[0,:] = x0 # Initial condition gives solution at first t + for i in range(len(t)-1): # Calculation loop Dt = t[i+1]-t[i] dxdt = array(f(x[i,:], t[i], *p)) - x[i+1,:] = x[i,:] + dxdt*Dt # Approximate solution at next value of x + # Approximate solution at next value of x + x[i+1,:] = x[i,:] + dxdt*Dt if verbose: - print('Numerical integration of ODE using explicit first-order method (Euler / Runge-Kutta) was successful.') + print('Numerical integration of ODE using explicit ' + + 'first-order method (Euler / Runge-Kutta) was successful.') return x def e2(f, x0, t, *p, verbose=False): @@ -143,19 +152,22 @@ def e2(f, x0, t, *p, verbose=False): :type x0: list :param t: time :type t: list - :param `*p`: parameters of the function (thickness, diameter, ...) + :param `*p`: parameters of the function (thickness, diameter, + ...) :param verbose: print information (default = False) :type verbose: bool """ - x = zeros((len(t), len(x0))) # Preallocate array - x[0,:] = x0 # Initial condition gives solution at first t - for i in range(len(t)-1): # Calculation loop + x = zeros((len(t), len(x0))) # Preallocate array + x[0,:] = x0 # Initial condition gives solution at first t + for i in range(len(t)-1): # Calculation loop Dt = t[i+1]-t[i] k_1 = array(f(x[i,:], t[i], *p)) k_2 = array(f(x[i,:]+0.5*Dt*k_1, t[i]+0.5*Dt, *p)) - x[i+1,:] = x[i,:] + k_2*Dt # Approximate solution at next value of x + # Approximate solution at next value of x + x[i+1,:] = x[i,:] + k_2*Dt if verbose: - print('Numerical integration of ODE using explicit 2th-order method (Runge-Kutta) was successful.') + print('Numerical integration of ODE using explicit ' + + '2th-order method (Runge-Kutta) was successful.') return x def e4(f, x0, t, *p, verbose=False): @@ -167,7 +179,8 @@ def e4(f, x0, t, *p, verbose=False): :type x0: list :param t: time :type t: list - :param `*p`: parameters of the function (thickness, diameter, ...) + :param `*p`: parameters of the function (thickness, diameter, + ...) :param verbose: print information (default = False) :type verbose: bool """ @@ -179,70 +192,64 @@ def e4(f, x0, t, *p, verbose=False): k_2 = array(f(x[i,:]+0.5*Dt*k_1, t[i]+0.5*Dt, *p)) k_3 = array(f(x[i,:]+0.5*Dt*k_2, t[i]+0.5*Dt, *p)) k_4 = array(f(x[i,:]+k_3*Dt, t[i]+Dt, *p)) - x[i+1,:] = x[i,:] + 1./6*(k_1+2*k_2+2*k_3+k_4)*Dt # Approximate solution at next value of x + # Approximate solution at next value of x + x[i+1,:] = x[i,:] + 1./6*(k_1+2*k_2+2*k_3+k_4)*Dt if verbose: - print('Numerical integration of ODE using explicit 4th-order method (Runge-Kutta) was successful.') + print('Numerical integration of ODE using explicit ' + + '4th-order method (Runge-Kutta) was successful.') return x -def dxdt_Dt(f, x, t, Dt, *p): - r""" - :param f: :math:`f = \dot{x}` - :type f: function - :param Dt: :math:`\Delta{t}` +def fpi(f, xi, ti, ti1, *p, max_iterations=1000, tol=1e-9, + verbose=False): + r"""Fixed-point iteration. - :returns: :math:`\Delta x = \dot{x} \Delta t` - """ - return array(dxdt(x, t, *p)) * Dt - -def fixed_point_iteration(f, xi, t, max_iterations=1000, tol=1e-9, verbose=False): - r""" - :param f: the function to iterate :math:`f = \Delta{x}(t)` + :param f: the function to iterate :math:`f = \dot{x}(x,t)` :type f: function :param xi: initial condition :math:`x_i` :type xi: list - :param t: time :math:`t` - :type t: float - :param `*p`: parameters of the function (thickness, diameter, ...) + :param ti: time :math:`t_i` + :type ti: float + :param ti1: time :math:`t_{i+1}` + :type ti1: float + :param `*p`: parameters of the function (thickness, diameter, + ...) :param max_iterations: maximum number of iterations :type max_iterations: int - :param tol: tolerance against residuum (default = 1e-9) + :param tol: tolerance against residuum :math:`\varepsilon` + (default = 1e-9) :type tol: float :param verbose: print information (default = False) :type verbose: bool - :returns: :math:`x_{i+1}` + :returns: :math:`x_{i}` .. math :: - x_{i+1} = x_i + \Delta x - - .. seealso:: - :meth:`dxdt_Dt` for :math:`\Delta x` + x_{i,j=0} = x_{i} + + .. math :: + x_{i,j+1} = x_i + \dot{x}(x_{i,j}, t_{i+1})\cdot(t_{i+1}-t_i) + + .. math :: + \text{residuum} = \frac{\lVert x_{i,j+1}-x_{i,j}\rVert} + {\lVert x_{i,j+1} \rVert} < \varepsilon + + .. math :: + x_{i} = x_{i,j=\text{end}} """ - xi = x0 - for j in range(max_iterations): # Fixed-point iteration - Dx = array(f(xi, t, *p)) - xi1 = x0 + Dx # Approximate solution at next value of x - residuum = norm(xi1-xi)/norm(xi1) - xi = xi1 + xij = xi + for j in range(max_iterations): + dxdt = array(f(xij, ti1, *p)) + # Approximate solution at next value of x + xij1 = xi + dxdt * (ti1-ti) + residuum = norm(xij1-xij)/norm(xij1) + xij = xij1 if residuum < tol: break iterations = j+1 # number beginning with 1 therefore + 1 - return xi, iterations + return xij, iterations -def i1n(f, x0, t, *p, max_iterations=1000, tol=1e-9, verbose=False): - iterations = zeros((len(t), 1)) - x = zeros((len(t), len(x0))) # Preallocate array - x[0,:] = x0 # Initial condition gives solution at first t - for i in range(len(t)-1): - Dt = t[i+1]-t[i] - xi = x[i,:] - Dx = dxdt_Dt(f, xi, t[i+1], Dt, *p) - x[i+1,:], iterations[i] = fixed_point_iteration(Dx, xi, t, max_iterations, tol, verbose) - if verbose: - print('Numerical integration of ODE using implicite first-order method (Euler) was successful.') - return x, iterations - -def i1(f, x0, t, *p, max_iterations=1000, tol=1e-9, verbose=False): +def i1(f, x0, t, *p, max_iterations=1000, tol=1e-9, + verbose=False): r"""Implicite first-order method / backward Euler method. :param f: the function to solve @@ -251,7 +258,8 @@ def i1(f, x0, t, *p, max_iterations=1000, tol=1e-9, verbose=False): :type x0: list :param t: time :type t: list - :param `*p`: parameters of the function (thickness, diameter, ...) + :param `*p`: parameters of the function (thickness, diameter, + ...) :param max_iterations: maximum number of iterations :type max_iterations: int :param tol: tolerance against residuum (default = 1e-9) @@ -263,26 +271,24 @@ def i1(f, x0, t, *p, max_iterations=1000, tol=1e-9, verbose=False): """ iterations = zeros((len(t), 1)) x = zeros((len(t), len(x0))) # Preallocate array - x[0,:] = x0 # Initial condition gives solution at first t + x[0,:] = x0 # Initial condition gives solution at first t + # x(i+1) = x(i) + f(x(i+1), t(i+1)), exact value of + # f(x(i+1), t(i+1)) is not available therefore using + # Newton-Raphson method for i in range(len(t)-1): Dt = t[i+1]-t[i] xi = x[i,:] - # x(i+1) = x(i) + f(x(i+1), t(i+1)), exact value of f(x(i+1), t(i+1)) is not - # available therefor using Newton-Raphson method - for j in range(max_iterations): # Fixed-point iteration - dxdt = array(f(xi, t[i+1], *p)) - xi1 = x[i,:] + dxdt*Dt # Approximate solution at next value of x - residuum = norm(xi1-xi)/norm(xi1) - xi = xi1 - if residuum < tol: - break - iterations[i] = j+1 + xi, iteration = fpi(f, xi, t[i], t[i+1], *p, max_iterations, + tol, verbose) x[i+1,:] = xi + iterations[i] = iteration if verbose: - print('Numerical integration of ODE using implicite first-order method (Euler) was successful.') + print('Numerical integration of ODE using implicite ' + + 'first-order method (Euler) was successful.') return x, iterations -def newmark_newtonraphson(f, x0, xp0, xpp0, t, *p, gamma=.5, beta=.25, max_iterations=1000, tol=1e-9, verbose=False): +def newmark_newtonraphson(f, x0, xp0, xpp0, t, *p, gamma=.5, + beta=.25, max_iterations=1000, tol=1e-9, verbose=False): r"""Newmark method. :param f: the function to solve @@ -295,7 +301,8 @@ def newmark_newtonraphson(f, x0, xp0, xpp0, t, *p, gamma=.5, beta=.25, max_itera :type xpp0: list :param t: time :type t: list - :param `*p`: parameters of the function (thickness, diameter, ...) + :param `*p`: parameters of the function (thickness, diameter, + ...) :param gamma: newmark parameter for velocity (default = 0.5) :type gamma: float :param beta: newmark parameter for displacement (default = 0.25) @@ -311,9 +318,9 @@ def newmark_newtonraphson(f, x0, xp0, xpp0, t, *p, gamma=.5, beta=.25, max_itera x = zeros((len(t), len(x0))) # Preallocate array xp = zeros((len(t), len(xp0))) # Preallocate array xpp = zeros((len(t), len(xpp0))) # Preallocate array - x[0,:] = x0 # Initial condition gives solution at first t - xp[0,:] = xp0 # Initial condition gives solution at first t - xpp[0,:] = xpp0 # Initial condition gives solution at first t + x[0,:] = x0 # Initial condition gives solution at first t + xp[0,:] = xp0 # Initial condition gives solution at first t + xpp[0,:] = xpp0 # Initial condition gives solution at first t for i in range(len(t)-1): Dt = t[i+1]-t[i] @@ -326,7 +333,8 @@ def newmark_newtonraphson(f, x0, xp0, xpp0, t, *p, gamma=.5, beta=.25, max_itera j = 0 for j in range(max_iterations): # Fixed-point iteration #dxdt = array(f(t[i+1], x1, p)) - #x11 = x[i,:] + dxdt*Dt # Approximate solution at next value of x + # Approximate solution at next value of x + #x11 = x[i,:] + dxdt*Dt N, dN, dNp, dNpp = f(x1.reshape(-1,).tolist(), xp1.reshape(-1,).tolist(), xpp1.reshape(-1,).tolist(), @@ -335,7 +343,8 @@ def newmark_newtonraphson(f, x0, xp0, xpp0, t, *p, gamma=.5, beta=.25, max_itera print('divergiert') break - xpp11 = xpp1 - dot(inv(dNpp), (N + dot(dN, (x1-xi)) + dot(dNp, (xp1-xpi)))) + xpp11 = xpp1 - dot(inv(dNpp), (N + dot(dN, (x1-xi)) + \ + dot(dNp, (xp1-xpi)))) xp1 = xpi + Dt*( (1-gamma)*xppi + gamma*xpp11 ) x1 = xi + Dt*xpi + Dt**2*( (.5-beta)*xppi + beta*xpp11 ) @@ -349,11 +358,13 @@ def newmark_newtonraphson(f, x0, xp0, xpp0, t, *p, gamma=.5, beta=.25, max_itera xp[i+1,:] = xp1.reshape(-1,).tolist() x[i+1,:] = x1.reshape(-1,).tolist() if verbose: - print('Numerical integration of ODE using explicite newmark method was successful.') + print('Numerical integration of ODE using explicite ' + + 'newmark method was successful.') return x, xp, xpp, iterations # x = concatenate((x, xp, xpp), axis=1) -def newmark_newtonraphson_rdk(fnm, x0, xp0, xpp0, t, *p, gamma=.5, beta=.25, maxIterations=1000, tol=1e-9, verbose=False): +def newmark_newtonraphson_rdk(fnm, x0, xp0, xpp0, t, *p, gamma=.5, + beta=.25, max_iterations=1000, tol=1e-9, verbose=False): r"""Newmark method. :param f: the function to solve @@ -366,7 +377,8 @@ def newmark_newtonraphson_rdk(fnm, x0, xp0, xpp0, t, *p, gamma=.5, beta=.25, max :type xpp0: list :param t: time :type t: list - :param `*p`: parameters of the function (thickness, diameter, ...) + :param `*p`: parameters of the function (thickness, diameter, + ...) :param gamma: newmark parameter for velocity (default = 0.5) :type gamma: float :param beta: newmark parameter for displacement (default = 0.25) @@ -382,13 +394,14 @@ def newmark_newtonraphson_rdk(fnm, x0, xp0, xpp0, t, *p, gamma=.5, beta=.25, max x = zeros((len(t), len(x0))) # Preallocate array xp = zeros((len(t), len(xp0))) # Preallocate array xpp = zeros((len(t), len(xpp0))) # Preallocate array - x[0,:] = x0 # Initial condition gives solution at first t - xp[0,:] = xp0 # Initial condition gives solution at first t - xpp[0,:] = xpp0 # Initial condition gives solution at first t + x[0,:] = x0 # Initial condition gives solution at first t + xp[0,:] = xp0 # Initial condition gives solution at first t + xpp[0,:] = xpp0 # Initial condition gives solution at first t for i in range(len(t)-1): Dt = t[i+1]-t[i] - rm, rmx, rmxpp, rd, rdx, rdxp, rk, rkx, f = fnm(x[i,:], xp[i,:], xpp[i,:], t[i], *p) + rm, rmx, rmxpp, rd, rdx, rdxp, rk, rkx, f = fnm(x[i,:], + xp[i,:], xpp[i,:], t[i], *p) xi = x[i,:].reshape(3,1) xpi = xp[i,:].reshape(3,1) @@ -397,13 +410,16 @@ def newmark_newtonraphson_rdk(fnm, x0, xp0, xpp0, t, *p, gamma=.5, beta=.25, max xp1 = xpi xpp1 = xppi j = 0 - for j in range(maxIterations): # Fixed-point iteration + for j in range(max_iterations): # Fixed-point iteration #dxdt = array(f(t[i+1], x1, p)) - #x11 = x[i,:] + dxdt*Dt # Approximate solution at next value of x + # Approximate solution at next value of x + #x11 = x[i,:] + dxdt*Dt r = (rmx+rdx+rkx)*Dt**2./4 + rdxp*Dt/2 + rmxpp - rp = f - (rm + dot(rmx, (Dt*xpi+Dt**2./4*xppi)) - dot(rmxpp, xppi) + \ - rd + dot(rdx, (Dt*xpi+Dt**2./4*xppi)) + dot(rdxp, Dt/2*xppi) + \ + rp = f - (rm + dot(rmx, (Dt*xpi+Dt**2./4*xppi)) - \ + dot(rmxpp, xppi) + \ + rd + dot(rdx, (Dt*xpi+Dt**2./4*xppi)) + \ + dot(rdxp, Dt/2*xppi) + \ rk + dot(rkx, (Dt*xpi+Dt**2./4*xppi)) ) xpp11 = dot(inv(r), rp) xp1 = xpi + Dt*( (1-gamma)*xppi + gamma*xpp11 ) @@ -419,6 +435,7 @@ def newmark_newtonraphson_rdk(fnm, x0, xp0, xpp0, t, *p, gamma=.5, beta=.25, max xp[i+1,:] = xp1.reshape(-1,).tolist() x[i+1,:] = x1.reshape(-1,).tolist() if verbose: - print('Numerical integration of ODE using explicite newmark method was successful.') + print('Numerical integration of ODE using explicite ' + + 'newmark method was successful.') return x, xp, xpp, iterations # x = concatenate((x, xp, xpp), axis=1) diff --git a/src/time_of_day.py b/src/time_of_day.py new file mode 100644 index 0000000..53c35ae --- /dev/null +++ b/src/time_of_day.py @@ -0,0 +1,156 @@ +#!/usr/bin/env python +# -*- coding: utf-8 -*- +"""Calculate time. + +:Date: 2019-06-01 + +.. module:: time_of_day + :platform: *nix, Windows + :synopsis: Calculate time. + +.. moduleauthor:: Daniel Weschke +""" +from __future__ import division, print_function, unicode_literals +from time import struct_time, mktime + + +def in_seconds(time): + """If time is `time.struct_time` convert to float seconds. + + :param time: the time in seconds + :type time: float or `time.struct_time` + + :returns: the time in seconds + :rtype: float + """ + if isinstance(time, struct_time): + time = mktime(time) + return time + +def seconds(time): + """The seconds of the time. + + :param time: the time in seconds + :type time: float or `time.struct_time` + + :returns: seconds, range [0, 60] + :rtype: float + """ + return in_seconds(time)%60 + +def seconds_norm(time): + """The seconds normalized to 60 seconds. + + :param time: the time in seconds + :type time: float or `time.struct_time` + + :returns: the normalized seconds, range [0, 1] + :rtype: float + """ + return seconds(time)/60 + +def minutes(time): + """The minutes of the time. + + :param time: the time in seconds + :type time: float or `time.struct_time` + + :returns: minutes, range [0, 60] + :rtype: float + """ + return in_seconds(time)/60%60 + +def minutes_norm(time): + """The minutes normalized to 60 minutes. + + :param time: the time in seconds + :type time: float or `time.struct_time` + + :returns: the normalized minutes, range [0, 1] + :rtype: float + """ + return minutes(time)/60 + +def hours(time): + """The hours of the time. + + :param time: the time in seconds + :type time: float or `time.struct_time` + + :returns: hours, range [0, 24] + :rtype: float + """ + return in_seconds(time)/60/60%24 + +def hours_norm(time): + """The hours normalized to 24 hours. + + :param time: the time in seconds + :type time: float or `time.struct_time` + + :returns: the normalized hours, range [0, 1] + :rtype: float + """ + return hours(time)/24 + +def days(time): + """The days of the time (year). + + :param time: the time in seconds + :type time: float or `time.struct_time` + + :returns: hours, range [0, 365.2425] + :rtype: float + """ + return in_seconds(time)/60/60/24%365.2425 + +def days_norm(time): + """The days normalized to 365.2425 (Gregorian, on average) days. + + :param time: the time in seconds + :type time: float or `time.struct_time` + + :returns: the normalized days, range [0, 1] + :rtype: float + """ + return days(time)/365.2425 + +def transform(time_norm, length, offset=0): + """Transform normalized time value to new length. + + :param position_norm: the normalized time value to transform + :type position_norm: float + :param length: the transformation + :type length: float + :param offset: the offset (default = 0) + :type offset: float + + :returns: the transformation value + :rtype: float + """ + return time_norm*length + offset + +if __name__ == "__main__": + from time import time, gmtime, localtime + # time in seconds + t = time() + min = minutes(t) + h = hours(t) + min_norm = minutes_norm(t) + h_norm = hours_norm(t) + + print('min ', min) + print('h ', h) + print('min_norm ', min_norm) + print('h_norm ', h_norm) + + x_len = 30 + x_offset = -8 + x_pos = transform(min_norm, x_len, x_offset) + print('m[-8,22] ', x_pos) + + y_len = 20 + y_offset = -10 + y_pos = transform(h_norm, y_len, y_offset) + print('h[-10,10]', y_pos) + diff --git a/tests/test_time_of_day.py b/tests/test_time_of_day.py new file mode 100644 index 0000000..12a2d2e --- /dev/null +++ b/tests/test_time_of_day.py @@ -0,0 +1,95 @@ +"""Test of fit module. + +:Date: 2019-06-01 + +.. module:: test_fit + :platform: *nix, Windows + :synopsis: Test of fit module. + +.. moduleauthor:: Daniel Weschke +""" +from __future__ import division, print_function, unicode_literals +import unittest +from time import mktime + +import os +import sys +sys.path.insert(0, os.path.abspath('../src')) +from time_of_day import (in_seconds, seconds, seconds_norm, + minutes, minutes_norm, hours, hours_norm, days, days_norm, + transform) + +# tm_year, tm_mon, tm_mday, tm_hour, tm_min, tm_sec, +# tm_wday, tm_yday, tm_isdst +T = mktime((2000, 11, 30, 0, 10, 20, 3, 335, -1)) + + +class TestTimeOfDay(unittest.TestCase): + + def test_in_seconds(self): + """test in_seconds""" + t = in_seconds(T) + self.assertEqual(t, 975539420.0) + + def test_seconds(self): + """test seconds""" + t = seconds(T) + self.assertEqual(t, 20.0) + + def test_seconds_norm(self): + """test seconds_norm""" + t = seconds_norm(T) + self.assertEqual(t, 20.0/60) # 0.3333333333333333 + + def test_minutes(self): + """test minutes""" + t = minutes(T) + self.assertEqual(t, T/60%60) # 10.333333333954215 + + def test_minutes_norm(self): + """test minutes_norm""" + t = minutes_norm(T) + self.assertEqual(t, T/60%60/60) # 0.17222222223257025 + + def test_hours(self): + """test hours""" + t = hours(T) + self.assertEqual(t, T/60/60%24) # 23.172222222259734 + + def test_hours_norm(self): + """test hours_norm""" + t = hours_norm(T) + self.assertEqual(t, T/60/60%24/24) # 0.9655092592608222 + + def test_days(self): + """test days""" + t = days(T) + self.assertEqual(t, T/60/60/24%365.2425) # 333.69050925926 + + def test_days_norm(self): + """test days_norm""" + t = days_norm(T) + # 0.9136135834664915 + self.assertEqual(t, T/60/60/24%365.2425/365.2425) + + def test_transform_minutes(self): + """test transform minutes""" + min_norm = minutes_norm(T) + x_len = 30 + x_offset = -8 + x_pos = transform(min_norm, x_len, x_offset) + # -2.8333333330228925 + self.assertEqual(x_pos, min_norm*x_len + x_offset) + + def test_transform_hours(self): + """test transform hours""" + h_norm = hours_norm(T) + y_len = 20 + y_offset = -10 + y_pos = transform(h_norm, y_len, y_offset) + # 9.310185185216444 + self.assertEqual(y_pos, h_norm*y_len + y_offset) + + +if __name__ == '__main__': + unittest.main()