========================================================== How to represent the DAE (Differential-Algebraic Equation) ========================================================== To compute consistent initial values of a DAE and to determine the index (and the corresponding condition number) you have to realize the DAE as an example class. Such a class object has the following attributes: * example_name - string delivering a description of the DAE * function_u - *True* or *False* depending on whether an additional function *u* is used or not. If the attribute *function_u* is absent, this corresponds to *False*. The function *u* can be used to fix initial values for some components or/and to realize (nonlinear) relations between components. See an example below. * type_of_DAE - string: 'proper', 'standard' or 'linear' :param proper: f and d are given for f(d'(x,t),x,t) :param standard: f is given for f(x',x,t) :param linear: f is given as A(t)x'+B(t)x - q(t) where A(t) and B(t) are realized as the functions eval_A and eval_B. * svdtol - double, which determines the zero level in rank determination. It has to be specified only in sensitive cases. The method uses the package **algopy** for algorithmic differentiation (AD), therefore you have to apply the methods from algopy to code the DAEs, i.e., all used functions like sin, cos, etc. and the vector/matrices like array, zeros, etc. have to be imported from algopy, see also examples below. ----------------------------------------------- An example written in all accepted formulations ----------------------------------------------- We give an example in standard formulation .. literalinclude:: example_Kronecker_standard.py formulated as proper DAE .. literalinclude:: example_Kronecker_proper.py and using the coefficient matrices for linear DAEs .. literalinclude:: example_Kronecker_linear.py ------------------------------------- A main program ------------------------------------- After coding the example we compute consistent initial values as well as the index of the DAE and its condition number at the computed point by calling the module *consistentValuesIndex* with the following input and output: :param: * **example** - object of example class * **TaylorOrder** - integer, order of the highest derivative of Taylor coefficients * **alpha** - reference vector as ndarray of size dim = dimension of the DAE * **x0** - ndarray containing initial guess of length dim*TaylorOrder * **t0** - timepoint * **augmentation** = False * if False: using the TaylorOrder * if True: augmentation of number of Taylor coefficients starting with 1 up to TaylorOrder+1 * **ftol** = 1e-8 - (parameter for python module *minimize*) * **disp** = False - (parameter for python module *minimize*) If *disp* = True also results of the iteration of the minimizer are displayed. * **maxiter** = 100 - (parameter for python module *minimize*): maximal number of iterations * **info_output** = 'ab' - a and b between 0 ... 3, higher level implies more information If a or b > 0 then two log-files are generated: *example_name*\ _array_\ *date*\ .log and *example_name*\ _consi_\ *date*\ .log \n *array* contains Jacobian matrices, projectors, ranks, SVD values etc.\n *consi* contains values related to the computation of consistent initial values and the solution in full accuracy (16 digits) :return: * **solution** as ndarray of length dim*(TaylorOrder+1) containing values for x(t0), x'(t0),x"(t0)/2!, ..., x :sup:`(TaylorOrder)`\ (t0)/(TaylorOrder)! **Caution**: Only x(t0), x'(t0),x"(t0)/2!, ..., x :sup:`(s)`\ (t0)/s! are consistent for s = TaylorOrder+1-index (see [1]) * **index** of the DAE in the computed consistent point * **cond** - condition number of the index computation * **computation_robust** - If True -> The relation TaylorOrder+1-index-s==0 is fulfilled We emphasize that the parameters **augmentation**, **ftol**, **disp**, **maxiter**, and **info_level** are optional, but **disp** = *True* and **info_level** > 0 should be set in order to obtain information. Here we give an example .. literalinclude:: Main_Kronecker_proper.py The terminal output is: .. literalinclude:: Result_MainKroneckerProper.txt This DAE in Kronecker canonical form has index 4 and is extensively discussed in [1]. In particular, it is explained there why only the values of the first s = TaylorOrder+1-index Taylor coefficients are consistent. --------------------------------------------------------------------- Pendulum example without/with *u* --------------------------------------------------------------------- We illustrate the use of the optional function *u* that prescribes additional requirements by means of a well-known nonlinear example. .. literalinclude:: example_pendulum_u_main.py The terminal output for function_u = False is: .. literalinclude:: ResultPendulumWithoutU.txt The terminal output for function_u = True is: .. literalinclude:: ResultPendulumWithU.txt ------------------------------------------------ Caraxis example without/with svdtol ------------------------------------------------ To determine the index of the DAE according to [2], Section 4, several rank considerations are made computing the corresponding singular values. In case the default tolerance does not lead to a robust rank-decision, the parameter svdtol should be set correspondingly in the example. We illustrate this with a well know example from `Test Set for IVP Solvers `_ .. literalinclude:: example_car_axis_.py If this example is run without setting svdtol, then the output reads .. literalinclude:: Result_example_car_axis.txt Note that, apart from this index-determination, we also check the s-fullness according to Corollaries 3 and 4 from [1]. If for the index mu and D=K+1 the relationship D-mu == s is not fulfilled, then some rank determinations may be critical. A clear hint that the rank-decision was not robust can be found if info_output[0] > 1 in the terminal output or in the related array-file (here e.g. car_axis_index3_array_*date*.log) .. literalinclude:: array_file_without_svdtol.txt Compare the latter singular value, which was considered to be nonzero (S(r) = 7.473449e-14), with the tolerance tol = 9.810314e-15. In contrast, if svdtol=1e-14, then a robust index characterization is obtained: .. literalinclude:: array_file_with_svdtol.txt In this case, the index 3 is correctly determined.