RK-Opt User Manual
Contents
1. 391 The upper and lower bounds on the unknowns These are chosen somewhat arbitrarily but usually aren t impor tant as long as they re not too restrictive 2 1 4 nonlinear constraints function con coneg nonlinear constraints x class s p objective poly coeff ind poly coeff val k Impose nonlinear constraints if objective ssp both order conditions and absolute monotonicity conditions if objective acc order conditions The input arguments are e vector of the decision variables See unpack rk m for details about the order in which they are stored class class of method to search erk explicit RK irk implicit RK dirk diagonally implicit RK sdirk singly diagonally implicit RK 25 35 2S 3S low storage formulations s number of stages p order of the RK scheme e objective objective function ssp maximize SSP coefficient minimize leading truncation error coefficient e poly ind index of the polynomial coefficients 5 for j gt p e poly val values of the polynomial coefficients 8 for j gt p tall tree elementary weights The outputs are con inequality constraints i e absolute monotonicity conditions if objective ssp or nothing if objective acc coneq order conditions plus stability function coefficients constraints tall tree elementary weights 8 C2. 5 Thus min s gt max p The equality holds for any order of accuracy because the number of linear order conditions that will be imposed to construct the GLM coefficients is p 2 2 2 radimpfast function rad radimpfast Compute the radius of absolute monotonicity of a rational function This function is outdated and needs to be fixed Uses van de Griend s algorithm assuming multiplicity one for all roots Uses high precision arithmetic 2 2 3 Rkp function R alpha beta Rkp k p Find the optimal SSP k step explicit LMM with order of accuracy p Inputs e k of steps 2 2 am_radius opt 13 RK Opt User Manual Release 0 2 order of accuracy Outputs a D the coefficients of the method Requires MATLAB s optimization toolbox for the LP solver 2 2 4 Rkp_dw function R alpha beta tbeta Rkp dw k p Finds the optimal SSP k step explicit LMM with order of accuracy p allowing downwind operators Inputs of steps p order of accuracy Outputs alpha beta tbeta the coefficients of the method The method is given by un zm alj thetalj tF Un n where tF u is the negated downwind operator Depends on MATLAB s optimization toolbox for the LP solver 2 2 5 Rkp imp function R alpha beta Rkp imp k p Find the optimal SSP k step implicit LMM with order of accuracy p Inputs e k of steps p order of accuracy Outputs alpha beta the coeffic
3. Some general utilities for analyzing Runge Kutta methods Contents RKtools am radius internal stab explicit butcher L2 timestep poly optimal shuosher form plotstabreg rk plotbounds s lw plotstabreg func rk stabfun semispectrum 2 4 1 am radius function r am radius A b c eps rmax Evaluates the Radius of absolute monotonicity of a Runge Kutta method given the Butcher array For an m stage method A should be an m x m matrix and b should be a column vector of length m Accuracy can be changed by modifying the value of eps default 107 10 Methods with very large radii of a m gt 50 will require the default value of rmax to be increased 16 Chapter 2 Reference RK Opt User Manual Release 0 2 The radius of absolute monotonicity is the largest value of r such that K I rA 0 2 1 tem lt m41 2 2 where A FF 2 4 2 internal_stab_explicit_butcher function stability internal stab explicit butcher A b c spectrum one step dt p This function computes and plots both intermediate and one step internal stability vector of an explicit Runge Kutta scheme given its Butcher tableau Note that for an explicit Runge Kutta scheme the stability functions are polynomials in the complex variable z Construct the intermediate stability functions psi j where j is the index of the stage Note that for an explicit scheme the intermediate stability polynomial
4. associated to the first stage is always 1 i e psi 1 1 Therefore we just compute and plot the remaining s 1 intermediate stability polynomials plus the one step stability polynomial of the Runge Kuatta method 2 4 3 L2 timestep poly function L2 timestep poly sdisc p d eps tol Find the absolutely timestep for a given combination of linear spatial discretization and stability function Also optionally plots the region of absolute stability and the eigenvalues The timestep is determined to within accuracy eps default 10 4 The spectral stability condition is checked to within tol default 10 13 2 4 4 optimal shuosher form function v alpha beta optimal shuosher form A b c 2 4 5 plotstabreg rk plotbounds ls lw function plotstabreg rk plotbounds ls lw Plots the absolute stability region of a Runge Kutta method given the Butcher array 2 4 6 plotstabreg func function dummy plotstabreg func p q bounds ls lw plot the absolute stability region of a one step method given the stability function Inputs 2 4 RKtools 17 RK Opt User Manual Release 0 2 p coefficients of the numerator of the stability function q coefficients of the denominator of the stability function if q is omitted it is assumed that the function is a polynomial Remaining inputs are optional bounds bounds for region to compute and plot default 9 1 5 5 e ls line style default r lw line width defa
5. opt 1 3 2 Running the tests To run the tests do the following in MATLAB gt gt cd path to RK Opt general test gt gt runLests If everything is set up correctly this will run several tests and inform you that the tests passed At present the tests are not very extensive 1 4 Citing RK Opt Are you using RK Opt in research work to be published If so please include explicit mention of our work in your publication We suggest language such as this 33 To solve problem 17 we used RK Opt a package for optimization of Runge Kutta methods 1 2 with the following entry in your bibliography 1 RK Opt Software for the design of Runge Kutta methods version 0 2 DI Ketcheson M Parsani and AJ Ahmadia http numerics kaust edu sa RK opt April 2013 as well as one or more of those below Specifically if you use the RK coeff opt package to optimize SSP coefficients please reference 2 If you use the RK coeff opt package to develop low storage methods please reference 3 If you use the RK coeff opt package to optimize for accuracy and or enforce a given stability function please reference 4 If you use the RK coeff opt package to develop general linear methods please reference 7 If you use the polyopt package please reference 5 If you use the am rad opt package please reference 6 4 Chapter 1 RK opt RK Opt User Manual Release 0 2 2 Highly Efficient Strong Stability Preserving Runge Kutt
6. 3 1773 2010 emsrk1 2 Explicit multistep Runge Kutta methods imsrk1 2 Implicit multistep Runge Kutta methods dimsrk1 2 Diagonally implicit multistep Runge Kutta methods objective objective function ssp maximize SSP coefficient acc minimize leading truncation error coefficient Accuracy optimization is not currently supported for multistep RK methods poly coeff ind index of the polynomial coefficients to constrain 5 for j gt p j denotes the index of the stage The default value is an empty array Note that one should not include any indices 7 lt p since those are determined by the order conditions poly coeff val constrained values of the polynomial coefficients 8 for 7 gt tall tree elementary weights The default value is an empty array startvec vector of the initial guess random random approach smart smart ap proach alternatively the user can provide the startvec array By default startvec is initial ize with random numbers solveorderconditions if set to 1 solve the order conditions first before trying to optimize The default value is 0 np number of processor to use If np gt 1 the MATLAB global optimization toolbox Multistart is used The default value is 1 just one core max tries maximum number of fmincon function calls The default value is 10 writeToFile whether to write to a file If set to 1 write the RK coefficient
7. 35 2S 38 low storage formulations s number of stages p order of the RK scheme e objective objective function ssp maximize SSP coefficient minimize leading truncation error coefficient The meaning of the output arguments is as follow 2 1 RK coeff opt 9 RK Opt User Manual Release 0 2 itis a scalar containing the radius of absolute monotonicity if objective ssp or the value of the leading truncation error coefficient if objective acc g itis a vector and contains the gradient of the objective function respect to the unknowns It is an array with all zero elements except for the last component which is equal to one if objective ssp or itis an empty array if objective 2 1 10 rk opt function rk rk opt s p class objective varargin Find optimal RK and multistep RK methods The meaning of the arguments is as follows s number of stages k number of steps 1 for RK methods p order of the Runge Kutta RK scheme class class of method to search Available classes erk Explicit Runge Kutta methods irk Implicit Runge Kutta methods dirk Diagonally implicit Runge Kutta methods sdirk Singly diagonally implicit Runge Kutta methods 2S etc Low storage explicit methods see Ketcheson Runge Kutta methods with minimum storage implementations J Comput Phys 229 5 176
8. RK Opt User Manual Release 0 2 David Ketcheson Matteo Parsani and Aron Ahmadia April 25 2013 CONTENTS 1 RK opt 3 1 1 Automated design of Runge Kutta methods gt s o ERR REOR 3 1 2 installation sns 4 66 Sandy Ie BES due ee eee Pie ERE ar qbus 3 1 3 Testing your installavion csse Quo hoes hee web eee bie eb EE S 4 LA RECS md ded Se rg 4 1 5 JAppliCadtofS ione e 8 See ee MER ESR Eom EORR e Der ew E E EY 5 2 Reference 7 21 JRK eoefbopt recitere o be bea ee hehe S OR RR S REO hae 7 2 2 am Tadius opt gage m bord ud OV IR RE uS dee e 12 2 9 roy GRY A EUNTEM Be See onum Dep DET S eps 15 2A 2254 99 ode ed Uere ete Peed weed Ew Rute e dei x eur 16 3 Contributing 19 RK Opt User Manual Release 0 2 RK opt is a MATLAB package for designing Runge Kutta RK methods and stability polynomials Supported ob jective functions include the principal error norm and the SSP coefficient Supported constraints include stability polynomial coefficients low storage formulations and structural constraints explicit diagonally implicit etc RK opt uses MATLAB s optimization toolbox in particular fmincon and linprog MATLAB s global optimization toolbox function Multistart can be used to exploit the benefits of parallel search on multicore machines The RK Opt package consists of the following packages RK
9. a Methods with Low Storage Imple mentations David I Ketcheson SIAM Journal on Scientific Computing 30 4 2113 2136 2008 3 Runge Kutta methods with minimum storage implementations David I Ketcheson Journal of Computational Physics 229 5 1763 1773 2010 4 Optimized explicit Runge Kutta schemes for the spectral difference method applied to wave propagation problems Matteo Parsani David I Ketcheson and W Deconinck SIAM Journal on Scien tific Computing In press Pre print available at http arxiv org abs 1207 5830 2013 5 Optimal stability polynomials for numerical integration of initial value problems David I Ketcheson and Aron J Ahmadia Communications in Applied Mathematics and Computational Science 7 2 247 271 2012 6 Computation of optimal monotonicity preserving general linear methods David I Ketcheson Mathematics of Computation 78 267 1497 1513 2009 7 Strong Stability Preserving Two step Runge Kutta Methods DI Ketcheson S Gottlieb CB Mac donald SIAM Journal on Numerical Analysis 2011 49 6 2618 2011 Also please do let us know if you are using this software so we can add your work to our Applications section 1 5 Applications RK Opt has been used in research for the following publications 1 Optimal implicit strong stability preserving Runge Kutta methods DI Ketcheson CB Macdon ald S Gottlieb Applied Numerical Mathematics 2009 59 2 373 392 2 Highly Efficient Strong Stab
10. asis varargin Finds an optimally stable polynomial of degree s and order p for the spectrum lam in the interval h min h max to precision eps Optional arguments lam func A function used to generate the appropriate spectrum at each bisection step instead of using a fixed scaled spectrum Used for instance to find the longest rectangle of a fixed height see Figure 10 of the CAMCOS paper Examples To find negative real axis inclusion lam spectrum realaxis 500 s 10 2 h poly_coeff opt poly bisect lam s p chebyshev To reproduce figure 10 of ketcheson ahmadia lam func spectrum rectangle 100 kappa 10 h poly coeff opt poly bisect lam 20 1 chebyshev lam func lam func plotstabreg func poly coeff 1 2 3 polyopt 15 RK Opt User Manual Release 0 2 2 3 2 spectrum function lamda spectrum name N kappa beta Return N discretely sampled values from certain sets in the complex plane Acceptable values for name e realaxis 1 0 jmagaxis 1 i disk z z 1 1 rectangle x lt lt lt lt 0 Niegemann ellipse and Niegemann circle See Niegemann 2011 gap Spectrum with a gap see Ketcheson amp Ahmadia 2012 kappa and beta are used only if name rectangle 2 3 3 test polyopt function test suite test polyopt A set of verification tests for the polyopt suite 2 4 RKtools
11. ckage 2 1 RK coeff opt This subpackage contains routines for finding optimal Runge Kutta method coefficents given a prescribed order of accuracy number of stages and an objective function Constraints on the stability polynomial possibly obtained using polyopt or am radius opt can optionally be provided Contents RK coeff opt check order errcoeff constraints nonlinear constraints oc albrecht oc butcher ksrk order conditions rk obj rk opt set n shuosher2butcher test SSP unpack lsrk unpack msrk unpack rk writeField 2 1 1 check RK order function check RK order A b c Determines order of a RK method up to sixth order For an s stage method input A should be a s x s matrix b and c should be column vectors of length s RK Opt User Manual Release 0 2 2 1 2 errcoeff function D errcoeff A b c p Inputs A b c Butcher tableau p order of accuracy of the method Computes the norm of the vector of truncation error coefficients for the terms of order 1 elementary weight 1 density of the tree symmetry of the tree For now we just use Butcher s approach We could alternatively use Albrecht s 2 1 3 linear constraints function Aeq beg lb ub linear constraints s class objective k This sets up The linear constraints corresponding to the consistency conditions j b land
12. coeff opt Find optimal Runge Kutta method coefficients for a prescribed order of accuracy and num ber of stages The objective function can be chosen as either the SSP coefficient or the leading truncation error coefficient The method may be constrained to have a low storage implementa tion and or a prescribed stability polynomial Implicit and diagonally implicit methods can also be optimized am radius opt Find stability functions with optimal radius of absolute monotonicity This includes codes for optimizing stability functions of multistep multistage methods and even methods with downwind ing The optimization of rational functions is experimental polyopt Given a spectrum typically corresponding to a spatial semi discretization of a PDE find an optimal stability polynomial in terms of its coefficients These polynomial coefficients can then be used as input to RK coeff opt to find a corresponding Runge Kutta method RKtools Some general utilities for analyzing Runge Kutta methods RK opt has been developed by David Ketcheson primary developer and maintainer Matteo Parsani and Aron Ah madia It is released under a modified BSD License If you use RK Opt in published work please see Citing RK Opt CONTENTS 1 RK Opt User Manual Release 0 2 2 CONTENTS RK OPT 1 1 Automated design of Runge Kutta methods An s stage Runge Kutta method has roughly s coefficients roughly s 2 for
13. e methods and even methods with downwinding Generally the optimization problem is phrased as a sequence of linear programming feasibility problems For details see ketcheson2009 The optimization of rational functions is experimental 12 Chapter 2 Reference RK Opt User Manual Release 0 2 Contents am radius opt multi opt radimpfast Rkp Rkp dw Rkp imp Rkp imp dw Rskp 2 2 1 multi R opt function multi multi R opt k p class varargin This function is a script to run the routines Rskp Rkp dw Rkp imp or Rkp imp dw several times with different inputs in order to construct tables of optimal values like those that appear in ketcheson2009 different values of the input parameters i e k1 k2 KK T length k ith element of steps p pl p2 pP T P length p ith element order of accuracy and 51 2 sSS T S length s ith element of stages when optimal contractive k step s stage GLM are investigated The family of method to be considered is specified in the string class Note that in general 5 5 5 P Fixed the order of accuracy of the time integration scheme is usually in terested in understanding the behavior of the threshold factor R as a function of the number of stages There SEI 6603 fore for a fixed element of the array this function loops over the elements of the array
14. explicit methods which can be chosen so as to provide high accuracy stability or other properties Historically most interest in Runge Kutta methods has focused on methods using the minimum number of stages for a given order of accuracy However in the past few decades there has been increasing recognition that using extra stages can be worthwhile in order to improve other method properties Some areas where this is particularly useful are in the enhancement of linear and nonlinear stability properties the reduction of storage requirements and the design of embedded pairs Methods with dozens or even hundreds of stages are not unheard of At the same time most existing Runge Kutta methods have been designed by hand by researchers laboriously solving the order conditions When using extra stages the number of available parameters makes the selection of a near optimal choice by hand impossible and one resorts to computational optimization This leads to a different paradigm of numerical method design in which we use sophisticated numerical optimization algorithms to design sophisticated numerical integration algorithms It can be expected that this trend will accelerate in the future and perhaps one day simple manually constructed algorithms will be the exception RK Opt contains a set of tools for designing Runge Kutta methods in this paradigm It has been constructed mostly in the direct line of our research but we have made some effort t
15. hapter 2 Reference RK Opt User Manual Release 0 2 Two forms of the order conditions are implemented one based on Butcher s approach and one based on Albrecht s approach One or the other may lead to a more tractable optimization problem in some cases but this has not been explored carefully The Albrecht order conditions are implemented up to order 9 assuming a certain stage order while the Butcher order conditions are implemented up to order 9 but do not assume anything about the stage order Albrecht s approach is used by default 2 1 5 oc albrecht function coneq oc_albrecht A b c p Order conditions for SSP RK Methods This version is based on Albrecht s approach 2 1 6 oc_butcher function coneq oc_butcher A b c p Order conditions for RKMs This version is based on Butcher s approach Assumes gt 1 2 1 7 oc_ksrk function coneq oc ksrk A b D theta p Order conditions for multistep RK methods 2 1 8 order_conditions function tau order conditions x class s p Aeq This is just a small wrapper used when solveorderconditions 1 2 1 9 rk obj function r g rk_obj x class s p objective Objective function for RK optimization The meaning of the input arguments is as follow vector of the unknowns e class class of method to search erk explicit RK irk implicit dirk diagonally implicit RK sdirk singly diagonally implicit RK 2S
16. ients of the method Depends on MATLAB s optimization toolbox for the LP solver 2 2 6 Rkp imp dw function R alpha beta Rkp imp dw k p Finds the optimal k step implicit LMM with order of accuracy p allowing downwinding Inputs of steps p order of accuracy Outputs alpha beta tbeta the coefficients of the method Depends on MATLAB s optimization toolbox for the LP solver 2 2 7 Rskp function R gamma Rskp s k p Finds the optimal contractive k step s stage GLM with order of accuracy p for linear problems Inputs s of stages of steps p order of accuracy 14 Chapter 2 Reference RK Opt User Manual Release 0 2 Outputs threshold factor gamma coefficients of the polynomials for k 1 the resulting polynomial is 5 7 o 1 z R in general the resulting stability function is Fill in epends on MATLAB s optimization toolbox for the LP solver 2 3 polyopt Given a spectrum typically corresponding to a spatial semi discretization of a PDE finds an optimal stability poly nomial The polynomial coefficients can then be used as input to RK coeff opt to find a corresponding Runge Kutta method This is the implementation of the algorithm described in ketcheson ahmadia The code was written by Aron Ahmadia and David Ketcheson Contents polyopt opt poly bisect spectrum test polyopt 2 3 1 opt poly bisect function h poly_coeff opt poly bisect lam s p b
17. ility Preserving Runge Kutta Methods with Low Storage Imple mentations David I Ketcheson SIAM Journal on Scientific Computing 30 4 2113 2136 2008 3 Runge Kutta methods with minimum storage implementations David I Ketcheson Journal of Computational Physics 229 5 1763 1773 2010 4 Optimized explicit Runge Kutta schemes for the spectral difference method applied to wave propagation problems Matteo Parsani David I Ketcheson and W Deconinck SIAM Journal on Scien tific Computing In press Pre print available at http arxiv org abs 1207 5830 2013 5 Optimal stability polynomials for numerical integration of initial value problems David I Ketcheson and Aron J Ahmadia Communications in Applied Mathematics and Computational Science 7 2 247 271 2012 6 Computation of optimal monotonicity preserving general linear methods David I Ketcheson Mathematics of Computation 78 267 1497 1513 2009 7 Strong Stability Preserving Two step Runge Kutta Methods DI Ketcheson S Gottlieb CB Mac donald SIAM Journal on Numerical Analysis 2011 49 6 2618 2011 8 Step Sizes for Strong Stability Preservation with Downwind Biased Operators DI Ketcheson SIAM Journal on Numerical Analysis 2011 49 4 1649 1 5 Applications 5 RK Opt User Manual Release 0 2 6 Chapter 1 RK opt REFERENCE This section contains a compilation of the documentation of each function organized by subpa
18. o help others easily understand and use it We hope that you find it useful and that you will contribute any enhancements you may develop back to the project by sending us a pull request on Github 1 2 Installation This section describes how to obtain RK opt and test that it is working correctly 1 2 1 Dependencies MATLAB 7 X or greater MATLAB optimization toolbox 1 2 2 Optional dependencies MATLAB Global Optimization toolbox for multicore search e xUnit test framework for testing your installation available for http www mathworks com matlabcentral fileexchange 22846 xUnit requires MATLAB R2008a or later RK Opt User Manual Release 0 2 1 2 3 Obtaining RK opt Download https github com ketch RK opt downloads Orclone git clone https github com ketch RK opt git After unzipping cloning add the subdirectory RK opt RKtools to your MATLAB path 1 3 Testing your installation You can test your RK opt installation with xUnit 1 3 1 Installing xUnit The MATLAB xUnit test framework can be downloaded for free at http www mathworks com matlabcentral fileexchange 22846 after the link click on the button in the upper right corner An easy way to install xUnit without setting any environ ment variables is to add the following lines to your startup m file addpath path to matlab xunit matlab xunit addpath path to matlab xunit matlab xunit xunit addpath path to RK opt RK coeff
19. s to a file called ERK p s txt The default value is 1 10 Chapter 2 Reference RK Opt User Manual Release 0 2 e algorithm which algorithm to use in fmincon sqp interior point or active set By default sqp is used Note numerical experiments have shown that when the objective function is the mini mization of the leading truncation error coefficient the interior point algorithm performs much better than the sqp one display level of display of fmincon solver off iter notify or final The default value is notify e problem class class of problems for which the RK is designed linear or nonlinear problems This option changes the type of order conditions check i e linear or nonlinear order conditions controll The default value is nonlinear Note Only 5 p class and objective are required inputs the other arguments are parameter name value arguments to the input parser scheme Therefore they can be specified in any order Example gt gt rk rk opt 4 3 erk acc max tries 2 np l solveorderconditions 1l The fmincon options are set through the optimset that creates alters optimization options structure By default the following ad MaxFunEvals 1000000 TolCon z 1 e 13 TolFun 1 e 13 ToIX 1 e 13 Maxlter 10000 Diagnostics off DerivativeCheck off GradObj on if the objective is set equal to
20. ssp 2 1 11 set n function n set_n s class Set total number of decision variables 2 1 12 shuosher2butcher function A b c shuosher2butcher alpha beta Generate Butcher form of a Runge Kutta method given its Shu Osher or modified Shu Osher form For an m stage method a and should be matrices of dimension m 1 x m 2 1 RK coeff opt 11 RK Opt User Manual Release 0 2 2 1 13 test SSP function test suite test SSP A set of verification tests for the RK opt package Currently this tests SSP coefficient optimization and accuracy optimization but not under constraints on the stability polynomial 2 1 14 unpack Isrk function A b bhat c alpha beta gammal gamma2 gamma3 delta unpack lsrk X s class Extracts the coefficient arrays from the optimization vector This function also returns the low storage coefficients 2 1 15 unpack msrk function A Ahat b bhat D theta unpack msrk X s k class Extract the coefficient arrays from the optimization vector 2 1 16 unpack rk function A b c unpack rk X s class Extracts the coefficient arrays from the optimization vector The coefficients are tored in a single vector x as b c A is stored row by row 2 1 17 writeField function wf writeField writeFid name value 2 2 am radius opt Find stability functions with optimal radius of absolute monotonicity This includes codes for optimizing stability functions of multistep multistag
21. ult 2 2 4 7 rk stabfun function p q rk stabfun rk Outputs the stability function of a RK method The Butcher coefficients should be stored in rk A rk b p contains the coefficients of the numerator q contains the coefficients of the denominator Ly pz _ 2 A eb7 oe Ejga de I 2A 2 4 8 semispectrum function L semispectrum method order doplot nx cfl Plot spectra of various semi discretizations of the advection equation Current choices for method fourier Fourier spectral method chebyshev Chebyshev spectral method updiff Upwind difference operators linearized WENO e DG Discontinuous Galerkin method The value of order matters only for the updiff and DG methods and selects the order of accuracy in those cases 18 Chapter 2 Reference CONTRIBUTING If you wish to contribute we recommend that you fork the RK Opt GitHub repository implement your additions and issue a pull request You may also simply e mail a patch to us 19
Download Pdf Manuals
Related Search