{"id":786,"date":"2013-10-14T12:00:10","date_gmt":"2013-10-14T17:00:10","guid":{"rendered":"https:\/\/blogs.mathworks.com\/cleve\/?p=786"},"modified":"2016-05-24T08:50:57","modified_gmt":"2016-05-24T13:50:57","slug":"complex-step-differentiation","status":"publish","type":"post","link":"https:\/\/blogs.mathworks.com\/cleve\/2013\/10\/14\/complex-step-differentiation\/","title":{"rendered":"Complex Step Differentiation"},"content":{"rendered":"<div class=\"content\"><!--introduction--><p>Complex step differentiation is a technique that employs complex arithmetic to obtain the numerical value of the first derivative of a real valued analytic function of a real variable, avoiding the loss of precision inherent in traditional finite differences.<\/p><!--\/introduction--><h3>Contents<\/h3><div><ul><li><a href=\"#7731c7fc-dc0e-46d3-9cff-ff916bbdbd41\">Stimulation<\/a><\/li><li><a href=\"#8177a107-64f8-4899-be52-6dd13197eca0\">Lyness and Moler<\/a><\/li><li><a href=\"#095a4500-3e18-453b-91d1-857a9f169a73\">The Algorithm<\/a><\/li><li><a href=\"#487a1d22-7fe2-446e-a2a3-69ec283c2d97\">An Example<\/a><\/li><li><a href=\"#39747585-d159-4360-b8e4-98ac58087010\">Symbolic Toolbox<\/a><\/li><li><a href=\"#80706c14-58cd-4434-98fe-cf80a3d60271\">Error Plot<\/a><\/li><li><a href=\"#8a8fa6b1-56ff-4987-a9ee-ed91e29cf07b\">References<\/a><\/li><\/ul><\/div><h4>Stimulation<a name=\"7731c7fc-dc0e-46d3-9cff-ff916bbdbd41\"><\/a><\/h4><p>It took prodding from a couple of people to bring about this post.<\/p><p>Many months ago I met Ricardo (Pax) Paxson, head of our Computational Biology group, in the MathWorks cafeteria.<\/p><div><ul><li>Pax: Hey, Cleve, you ought to blog about that complex step algorithm you invented.<\/li><li>Me: What are you talking about?<\/li><li>Pax: You know, differentiation using complex arithmetic.  It's hot stuff. We're using it in SimBiology.<\/li><li>Me: I don't remember inventing any complex step differentiation algorithm.<\/li><li>Pax: Well, an old paper of yours is referenced.  Lots of people are using it today.<\/li><\/ul><\/div><p>Then last February I met Joaquim Martins, a Professor of Aerospace Engineering at the University of Michigan, during the SIAM CS&amp;E meeting. He gave me three papers he had written about the complex step method that he and his colleagues were using to compute derivatives of simulations that couple computational fluid dynamics with structural finite-element analysis for aircraft design applications<\/p><p>Joaquim's papers refer to a paper on numerical differentiation using complex arithmetic that James Lyness and I had published in the <i>SIAM Journal of Numerical Analysis<\/i> in 1967.  That's almost 50 years ago.<\/p><h4>Lyness and Moler<a name=\"8177a107-64f8-4899-be52-6dd13197eca0\"><\/a><\/h4><p>James Lyness was a buddy of mine.  We met at the ETH in Zurich when I was visiting there on my postdoc in 1965-66.  He was an Englishman a few years older than I, visiting the ETH at roughly the same time. After that visit, he spent the rest of his professional career at Argonne Laboratory.  He didn't work on LINPACK and EISPACK, but I would see him at Argonne when I visited there to work on those packages. The photo that I often show of the four LINPACK authors lounging on the back of Jack Dongarra's car was taken during a party at the Lyness home.<\/p><p>That year in Zurich was very productive research wise.  James had lots of interesting ideas and we published several papers.  One of them was \"Numerical Differentiation of Analytic Functions\".  It introduces, I guess for the first time, the idea of using complex arithmetic on functions of a real variable.  But the actual algorithms are pretty exotic. They are derived from Cauchy's Integral Theorem and involve things like Moebius numbers and the Euler Maclaurin summation formula.  If you asked me to describe them in any detail today, I wouldn't be able to.<\/p><p>The complex step algorithm that I'm about to describe is much simpler than anything described in that 1967 paper.  So, although we have some claim of precedence for the idea of using complex arithmetic on real functions, we certainly cannot claim to have invented today's algorithm.<\/p><p>The complex step algorithm was first described in a brief paper by William Squire and George Trapp in <i>SIAM Review<\/i> in 1998.  Martins, together with his colleagues and students, offers an extensive guide in the papers I refer to below.<\/p><h4>The Algorithm<a name=\"095a4500-3e18-453b-91d1-857a9f169a73\"><\/a><\/h4><p>We're concerned with an <i>analytic<\/i> function.  Mathematically, that means the function is infinitely differentiable and can be smoothly extended into the complex plane.  Computationally, it probably means that it is defined by a single \"one line\" formula, not a more extensive piece of code with <tt>if<\/tt> statements and <tt>for<\/tt> loops.<\/p><p>Let $F(z)$ be such a function, let $x_0$ be a point on the real axis, and let $h$ be a real parameter.  Expand $F(z)$ in a Taylor series off the real axis.<\/p><p>$$ F(x_0+ih) = F(x_0)+ihF'(x_0)-h^2F''(x_0)\/2!-ih^3F^{(3)}\/3!+... $$<\/p><p>Take the imaginary part of both sides and divide by $h$.<\/p><p>$$ F'(x_0) = Im(F(x_0+ih))\/h + O(h^2) $$<\/p><p>Simply evaluating the function $F$ at the imaginary argument $x_0+ih$, and dividing by $h$, gives an approximation to the value of the derivative, $F'(x_0)$, that is accurate to order $O(h^2)$.  We might as well choose $h=10^{-8}$.  Then the error in the approximation is about the same size as the roundoff error involved in storing a double precision floating point value of $F'(x_0)$.<\/p><p>Here, then, is the complex step differentiation algorithm in its entirety:<\/p><p>$$ h = 10^{-8} $$<\/p><p>$$ F'(x_0) \\approx Im(F(x_0+ih))\/h $$<\/p><h4>An Example<a name=\"487a1d22-7fe2-446e-a2a3-69ec283c2d97\"><\/a><\/h4><p>Here is one of my favorite examples.  I used it in the 1967 paper and it has become a kind of standard in this business ever since.<\/p><p>$$ F(x) = \\frac{\\mbox{e}^{x}}{(\\cos{x})^3+(\\sin{x})^3} $$<\/p><pre class=\"codeinput\">   F = @(x) exp(x).\/((cos(x)).^3 + (sin(x)).^3)\r\n<\/pre><pre class=\"codeoutput\">F = \r\n    @(x)exp(x).\/((cos(x)).^3+(sin(x)).^3)\r\n<\/pre><p>Here's a plot over a portion of the $x$-axis, with a dot at $\\pi\/4$.<\/p><pre class=\"codeinput\">   ezplot(F,[-pi\/4,pi\/2])\r\n   axis([-pi\/4,pi\/2,0,6])\r\n   set(gca,<span class=\"string\">'xtick'<\/span>,[-pi\/4,0,pi\/4,pi\/2])\r\n   line([pi\/4,pi\/4],[F(pi\/4),F(pi\/4)],<span class=\"string\">'marker'<\/span>,<span class=\"string\">'.'<\/span>,<span class=\"string\">'markersize'<\/span>,18)\r\n<\/pre><img decoding=\"async\" vspace=\"5\" hspace=\"5\" src=\"https:\/\/blogs.mathworks.com\/images\/cleve\/cstep_01.png\" alt=\"\"> <p>Let's find the slope at the dot, that's $F'(\\pi\/4)$.<\/p><p>The complex step approximation where I have not yet fixed the step size <tt>h<\/tt> is<\/p><pre class=\"codeinput\">   Fpc = @(x,h) imag(F(x+i*h))\/h;\r\n<\/pre><p>The centered finite difference approximation with step size <tt>h<\/tt> is<\/p><pre class=\"codeinput\">   Fpd = @(x,h) (F(x+h) - F(x-h))\/(2*h);\r\n<\/pre><p>Compare them for a sequence of decreasing <tt>h<\/tt>.<\/p><pre class=\"codeinput\">   format <span class=\"string\">long<\/span>\r\n   disp([<span class=\"string\">'           h             complex step     finite differerence'<\/span>])\r\n   <span class=\"keyword\">for<\/span> h = 10.^(-(1:16))\r\n      disp([h Fpc(pi\/4,h) Fpd(pi\/4,h)])\r\n   <span class=\"keyword\">end<\/span>\r\n<\/pre><pre class=\"codeoutput\">           h             complex step     finite differerence\r\n   0.100000000000000   3.144276040634560   3.061511866568119\r\n   0.010000000000000   3.102180075411270   3.101352937655877\r\n   0.001000000000000   3.101770529535847   3.101762258158169\r\n   0.000100000000000   3.101766435192940   3.101766352480162\r\n   0.000010000000000   3.101766394249620   3.101766393398542\r\n   0.000001000000000   3.101766393840188   3.101766393509564\r\n   0.000000100000000   3.101766393836091   3.101766394841832\r\n   0.000000010000000   3.101766393836053   3.101766443691645\r\n   0.000000001000000   3.101766393836052   3.101766399282724\r\n   0.000000000100000   3.101766393836052   3.101763290658255\r\n   0.000000000010000   3.101766393836053   3.101785495118747\r\n   0.000000000001000   3.101766393836053   3.101963130802687\r\n   0.000000000000100   3.101766393836052   3.097522238704187\r\n   0.000000000000010   3.101766393836053   3.108624468950438\r\n   0.000000000000001   3.101766393836052   3.108624468950438\r\n   0.000000000000000   3.101766393836053   8.881784197001252\r\n<\/pre><p>As predicted, complex step has obtained its full double precision result half way through the sequence of <tt>h's<\/tt>.  On the other hand, the finite difference method never achieves more than half double precision accuracy and completely breaks down by the time <tt>h<\/tt> reaches the last value.<\/p><h4>Symbolic Toolbox<a name=\"39747585-d159-4360-b8e4-98ac58087010\"><\/a><\/h4><p>Since we have access to the Symbolic Toolbox, we can get the exact answer. Redefine <tt>x<\/tt> and <tt>F<\/tt> to be symbolic.<\/p><pre class=\"codeinput\">   syms <span class=\"string\">x<\/span>\r\n   F = exp(x)\/((cos(x))^3 + (sin(x))^3)\r\n<\/pre><pre class=\"codeoutput\">F =\r\nexp(x)\/(cos(x)^3 + sin(x)^3)\r\n<\/pre><p>Differentiate and evaluate the derivative at $\\pi\/4$.<\/p><pre class=\"codeinput\">   Fps = simple(diff(F))\r\n   exact = subs(Fps,pi\/4)\r\n   flexact = double(exact)\r\n<\/pre><pre class=\"codeoutput\">Fps =\r\n(exp(x)*(cos(3*x) + sin(3*x)\/2 + (3*sin(x))\/2))\/(cos(x)^3 + sin(x)^3)^2\r\nexact =\r\n2^(1\/2)*exp(pi\/4)\r\nflexact =\r\n   3.101766393836051\r\n<\/pre><p>This verifies the result we obtained with the complex step.<\/p><h4>Error Plot<a name=\"80706c14-58cd-4434-98fe-cf80a3d60271\"><\/a><\/h4><p>Now that we now the exact answer, we can plot the error from computing the first derivative by complex step and finite differences.<\/p><pre class=\"codeinput\">   h = 10.^(-(1:16)');\r\n   errs = zeros(16,2);\r\n   <span class=\"keyword\">for<\/span> k = 1:16\r\n      errs(k,1) = abs(Fpc(pi\/4,h(k))-exact);\r\n      errs(k,2) = abs(Fpd(pi\/4,h(k))-exact);\r\n   <span class=\"keyword\">end<\/span>\r\n   loglog(h,errs)\r\n   axis([eps 1 eps 1])\r\n   set(gca,<span class=\"string\">'xdir'<\/span>,<span class=\"string\">'rev'<\/span>)\r\n   legend(<span class=\"string\">'complex step'<\/span>,<span class=\"string\">'finite difference'<\/span>,<span class=\"string\">'location'<\/span>,<span class=\"string\">'southwest'<\/span>)\r\n   xlabel(<span class=\"string\">'step size h'<\/span>)\r\n   ylabel(<span class=\"string\">'error'<\/span>)\r\n<\/pre><img decoding=\"async\" vspace=\"5\" hspace=\"5\" src=\"https:\/\/blogs.mathworks.com\/images\/cleve\/cstep_02.png\" alt=\"\"> <p>This confirms that complex step is accurate for small step sizes whereas the finite difference approach never achieves full accuracy.<\/p><h4>References<a name=\"8a8fa6b1-56ff-4987-a9ee-ed91e29cf07b\"><\/a><\/h4><p>Lyness, James N., and Moler, Cleve B., Numerical Differentiation of Analytic Functions, <i>SIAM J. Numerical Analysis 4<\/i>, 1967, pp. 202-210. <a href=\"http:\/\/epubs.siam.org\/doi\/abs\/10.1137\/0704019\">epubs.siam.org\/doi\/abs\/10.1137\/0704019<\/a><\/p><p>Squire, William, and Trapp, George, Using complex variables to estimate derivatives of real functions, <i>SIAM Review 40<\/i>, 1998, pp. 110-112. <a href=\"http:\/\/epubs.siam.org\/doi\/abs\/10.1137\/S003614459631241X\">epubs.siam.org\/doi\/abs\/10.1137\/S003614459631241X<\/a><\/p><p>Martins, Joaquim R. R. A., Sturdza, Peter, and Alonso, Juan J., The Complex-Step Derivative Approximation, <i>ACM Transactions on Mathematical Software 29<\/i>, 2003, pp. 245-262, <a href=\"\">Martins2003CSD.pdf<\/a><\/p><p>Martins, Joaquim R. R. A., A Guide to the Complex-Step Derivative Approximation, <a href=\"\">complex step guide<\/a><\/p><p>Kenway, Gaetan K.W., Kennedy, Graeme J., and Martins, Joaquim R. R. A., A scalable parallel approach for high-fidelity steady-state aeroelastic analysis and adjoint derivative computations, AIAA Journal (In Press)<\/p><script language=\"JavaScript\"> <!-- \r\n    function grabCode_5382da4c6b214fbf8340cdef7eab1dcf() {\r\n        \/\/ Remember the title so we can use it in the new page\r\n        title = document.title;\r\n\r\n        \/\/ Break up these strings so that their presence\r\n        \/\/ in the Javascript doesn't mess up the search for\r\n        \/\/ the MATLAB code.\r\n        t1='5382da4c6b214fbf8340cdef7eab1dcf ' + '##### ' + 'SOURCE BEGIN' + ' #####';\r\n        t2='##### ' + 'SOURCE END' + ' #####' + ' 5382da4c6b214fbf8340cdef7eab1dcf';\r\n    \r\n        b=document.getElementsByTagName('body')[0];\r\n        i1=b.innerHTML.indexOf(t1)+t1.length;\r\n        i2=b.innerHTML.indexOf(t2);\r\n \r\n        code_string = b.innerHTML.substring(i1, i2);\r\n        code_string = code_string.replace(\/REPLACE_WITH_DASH_DASH\/g,'--');\r\n\r\n        \/\/ Use \/x3C\/g instead of the less-than character to avoid errors \r\n        \/\/ in the XML parser.\r\n        \/\/ Use '\\x26#60;' instead of '<' so that the XML parser\r\n        \/\/ doesn't go ahead and substitute the less-than character. \r\n        code_string = code_string.replace(\/\\x3C\/g, '\\x26#60;');\r\n\r\n        copyright = 'Copyright 2013 The MathWorks, Inc.';\r\n\r\n        w = window.open();\r\n        d = w.document;\r\n        d.write('<pre>\\n');\r\n        d.write(code_string);\r\n\r\n        \/\/ Add copyright line at the bottom if specified.\r\n        if (copyright.length > 0) {\r\n            d.writeln('');\r\n            d.writeln('%%');\r\n            if (copyright.length > 0) {\r\n                d.writeln('% _' + copyright + '_');\r\n            }\r\n        }\r\n\r\n        d.write('<\/pre>\\n');\r\n\r\n        d.title = title + ' (MATLAB code)';\r\n        d.close();\r\n    }   \r\n     --> <\/script><p style=\"text-align: right; font-size: xx-small; font-weight:lighter;   font-style: italic; color: gray\"><br><a href=\"javascript:grabCode_5382da4c6b214fbf8340cdef7eab1dcf()\"><span style=\"font-size: x-small;        font-style: italic;\">Get \r\n      the MATLAB code <noscript>(requires JavaScript)<\/noscript><\/span><\/a><br><br>\r\n      Published with MATLAB&reg; R2013b<br><\/p><p class=\"footer\"><br>\r\n      Published with MATLAB&reg; R2013b<br><\/p><\/div><!--\r\n5382da4c6b214fbf8340cdef7eab1dcf ##### SOURCE BEGIN #####\r\n%% Complex Step Differentiation\r\n% Complex step differentiation is a technique that employs complex arithmetic\r\n% to obtain the numerical value of the first derivative of a real valued\r\n% analytic function of a real variable, avoiding the loss of precision\r\n% inherent in traditional finite differences.\r\n\r\n%% Stimulation\r\n% It took prodding from a couple of people to bring about this post.\r\n\r\n%%\r\n% Many months ago I met Ricardo (Pax) Paxson, head of our Computational\r\n% Biology group, in the MathWorks cafeteria.\r\n%\r\n% * Pax: Hey, Cleve, you ought to blog about that complex step algorithm\r\n% you invented. \r\n% * Me: What are you talking about?\r\n% * Pax: You know, differentiation using complex arithmetic.  It's hot stuff.\r\n% We're using it in SimBiology.\r\n% * Me: I don't remember inventing any complex step differentiation algorithm.\r\n% * Pax: Well, an old paper of yours is referenced.  Lots of people are\r\n% using it today.\r\n\r\n%%\r\n% Then last February I met Joaquim Martins, a Professor of Aerospace\r\n% Engineering at the University of Michigan, during the SIAM CS&E meeting.\r\n% He gave me three papers he had written about the complex step method that\r\n% he and his colleagues were using to compute derivatives of simulations\r\n% that couple computational fluid dynamics with structural finite-element\r\n% analysis for aircraft design applications\r\n\r\n%%\r\n% Joaquim's papers refer to a paper on numerical differentiation using\r\n% complex arithmetic that James Lyness and I had published in the\r\n% _SIAM Journal of Numerical Analysis_ in 1967.  That's almost 50 years ago.\r\n\r\n%% Lyness and Moler\r\n% James Lyness was a buddy of mine.  We met at the ETH in Zurich when I\r\n% was visiting there on my postdoc in 1965-66.  He was an Englishman a\r\n% few years older than I, visiting the ETH at roughly the same time.\r\n% After that visit, he spent the rest of his professional career at\r\n% Argonne Laboratory.  He didn't work on LINPACK and EISPACK, but I would\r\n% see him at Argonne when I visited there to work on those packages.\r\n% The photo that I often show of the four LINPACK authors lounging on the\r\n% back of Jack Dongarra's car was taken during a party at the Lyness home.\r\n\r\n%%\r\n% That year in Zurich was very productive research wise.  James had lots\r\n% of interesting ideas and we published several papers.  One of them was\r\n% \"Numerical Differentiation of Analytic Functions\".  It introduces, I guess\r\n% for the first time, the idea of using complex arithmetic on functions of\r\n% a real variable.  But the actual algorithms are pretty exotic.\r\n% They are derived from Cauchy's Integral Theorem and involve things like\r\n% Moebius numbers and the Euler Maclaurin summation formula.  If you asked me\r\n% to describe them in any detail today, I wouldn't be able to.\r\n\r\n%%\r\n% The complex step algorithm that I'm about to describe is much simpler\r\n% than anything described in that 1967 paper.  So, although we have some\r\n% claim of precedence for the idea of using complex arithmetic on real\r\n% functions, we certainly cannot claim to have invented today's algorithm.\r\n\r\n%%\r\n% The complex step algorithm was first described in a brief paper by\r\n% William Squire and George Trapp in _SIAM Review_ in 1998.  Martins,\r\n% together with his colleagues and students, offers an extensive guide in the\r\n% papers I refer to below.\r\n\r\n%% The Algorithm\r\n% We're concerned with an _analytic_ function.  Mathematically, that means\r\n% the function is infinitely differentiable and can be smoothly extended into\r\n% the complex plane.  Computationally, it probably means that it is defined\r\n% by a single \"one line\" formula, not a more extensive piece of code with\r\n% |if| statements and |for| loops.\r\n\r\n%%\r\n% Let $F(z)$ be such a function, let $x_0$ be a point on the real axis,\r\n% and let $h$ be a real parameter.  Expand $F(z)$ in a Taylor series off\r\n% the real axis.\r\n%\r\n% $$ F(x_0+ih) = F(x_0)+ihF'(x_0)-h^2F''(x_0)\/2!-ih^3F^{(3)}\/3!+... $$\r\n%\r\n% Take the imaginary part of both sides and divide by $h$.\r\n%\r\n% $$ F'(x_0) = Im(F(x_0+ih))\/h + O(h^2) $$\r\n%\r\n% Simply evaluating the function $F$ at the imaginary argument $x_0+ih$,\r\n% and dividing by $h$, gives an approximation to the value of the derivative,\r\n% $F'(x_0)$, that is accurate to order $O(h^2)$.  We might as well choose\r\n% $h=10^{-8}$.  Then the error in the approximation is about the same\r\n% size as the roundoff error involved in storing a double precision floating\r\n% point value of $F'(x_0)$.\r\n\r\n%%\r\n% Here, then, is the complex step differentiation algorithm in its entirety:\r\n%\r\n% $$ h = 10^{-8} $$\r\n%\r\n% $$ F'(x_0) \\approx Im(F(x_0+ih))\/h $$\r\n\r\n\r\n%% An Example\r\n% Here is one of my favorite examples.  I used it in the 1967 paper and it\r\n% has become a kind of standard in this business ever since.\r\n%\r\n% $$ F(x) = \\frac{\\mbox{e}^{x}}{(\\cos{x})^3+(\\sin{x})^3} $$\r\n%\r\n   F = @(x) exp(x).\/((cos(x)).^3 + (sin(x)).^3)\r\n\r\n%%\r\n% Here's a plot over a portion of the $x$-axis, with a dot at $\\pi\/4$.\r\n\r\n   ezplot(F,[-pi\/4,pi\/2])\r\n   axis([-pi\/4,pi\/2,0,6])\r\n   set(gca,'xtick',[-pi\/4,0,pi\/4,pi\/2])\r\n   line([pi\/4,pi\/4],[F(pi\/4),F(pi\/4)],'marker','.','markersize',18)\r\n\r\n%%\r\n% Let's find the slope at the dot, that's $F'(\\pi\/4)$.\r\n%\r\n%%\r\n% The complex step approximation where I have not yet fixed the step\r\n% size |h| is\r\n%\r\n   Fpc = @(x,h) imag(F(x+i*h))\/h;\r\n\r\n%%\r\n% The centered finite difference approximation with step size |h| is\r\n\r\n   Fpd = @(x,h) (F(x+h) - F(x-h))\/(2*h);\r\n\r\n%%\r\n% Compare them for a sequence of decreasing |h|.\r\n\r\n   format long\r\n   disp(['           h             complex step     finite differerence'])\r\n   for h = 10.^(-(1:16))\r\n      disp([h Fpc(pi\/4,h) Fpd(pi\/4,h)])\r\n   end\r\n\r\n%%\r\n% As predicted, complex step has obtained its full double precision result\r\n% half way through the sequence of |h's|.  On the other hand, the finite\r\n% difference method never achieves more than half double precision accuracy\r\n% and completely breaks down by the time |h| reaches the last value.\r\n\r\n%% Symbolic Toolbox\r\n% Since we have access to the Symbolic Toolbox, we can get the exact answer.\r\n% Redefine |x| and |F| to be symbolic.\r\n\r\n   syms x\r\n   F = exp(x)\/((cos(x))^3 + (sin(x))^3)\r\n\r\n%%\r\n% Differentiate and evaluate the derivative at $\\pi\/4$.\r\n\r\n   Fps = simple(diff(F))\r\n   exact = subs(Fps,pi\/4)\r\n   flexact = double(exact)\r\n\r\n%%\r\n% This verifies the result we obtained with the complex step.\r\n\r\n%% Error Plot\r\n% Now that we now the exact answer, we can plot the error from computing\r\n% the first derivative by complex step and finite differences.\r\n\r\n   h = 10.^(-(1:16)');\r\n   errs = zeros(16,2);  \r\n   for k = 1:16\r\n      errs(k,1) = abs(Fpc(pi\/4,h(k))-exact);\r\n      errs(k,2) = abs(Fpd(pi\/4,h(k))-exact);\r\n   end\r\n   loglog(h,errs)\r\n   axis([eps 1 eps 1])\r\n   set(gca,'xdir','rev')\r\n   legend('complex step','finite difference','location','southwest')\r\n   xlabel('step size h')\r\n   ylabel('error')\r\n   \r\n%%\r\n% This confirms that complex step is accurate for small step sizes whereas\r\n% the finite difference approach never achieves full accuracy.\r\n   \r\n\r\n%% References\r\n% Lyness, James N., and Moler, Cleve B.,\r\n% Numerical Differentiation of Analytic Functions,\r\n% _SIAM J. Numerical Analysis 4_, 1967, pp. 202-210. \r\n% <http:\/\/epubs.siam.org\/doi\/abs\/10.1137\/0704019\r\n% epubs.siam.org\/doi\/abs\/10.1137\/0704019>\r\n\r\n%%\r\n% Squire, William, and Trapp, George, \r\n% Using complex variables to estimate derivatives of real functions,\r\n% _SIAM Review 40_, 1998, pp. 110-112.\r\n% <http:\/\/epubs.siam.org\/doi\/abs\/10.1137\/S003614459631241X\r\n% epubs.siam.org\/doi\/abs\/10.1137\/S003614459631241X>\r\n\r\n%%\r\n% Martins, Joaquim R. R. A., Sturdza, Peter, and Alonso, Juan J.,\r\n% The Complex-Step Derivative Approximation,\r\n% _ACM Transactions on Mathematical Software 29_, 2003, pp. 245-262,\r\n% <\r\n% Martins2003CSD.pdf>\r\n\r\n%%\r\n% Martins, Joaquim R. R. A.,\r\n% A Guide to the Complex-Step Derivative Approximation,\r\n% <\r\n% complex step guide>\r\n\r\n%%\r\n% Kenway, Gaetan K.W., Kennedy, Graeme J., and Martins, Joaquim R. R. A.,\r\n% A scalable parallel approach for high-fidelity steady-state aeroelastic\r\n% analysis and adjoint derivative computations, \r\n% AIAA Journal (In Press), \r\n% <http:\/\/mdolab.engin.umich.edu\/sites\/default\/files\/kenway2013a_0.pdf\r\n% kenway2013a.pdf>\r\n\r\n\r\n##### SOURCE END ##### 5382da4c6b214fbf8340cdef7eab1dcf\r\n-->","protected":false},"excerpt":{"rendered":"<div class=\"overview-image\"><img decoding=\"async\"  class=\"img-responsive\" src=\"https:\/\/blogs.mathworks.com\/images\/cleve\/cstep_01.png\" onError=\"this.style.display ='none';\" \/><\/div><!--introduction--><p>Complex step differentiation is a technique that employs complex arithmetic to obtain the numerical value of the first derivative of a real valued analytic function of a real variable, avoiding the loss of precision inherent in traditional finite differences.... <a class=\"read-more\" href=\"https:\/\/blogs.mathworks.com\/cleve\/2013\/10\/14\/complex-step-differentiation\/\">read more >><\/a><\/p>","protected":false},"author":78,"featured_media":0,"comment_status":"open","ping_status":"closed","sticky":false,"template":"","format":"standard","meta":[],"categories":[11,12,4,16],"tags":[],"_links":{"self":[{"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/posts\/786"}],"collection":[{"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/posts"}],"about":[{"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/types\/post"}],"author":[{"embeddable":true,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/users\/78"}],"replies":[{"embeddable":true,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/comments?post=786"}],"version-history":[{"count":10,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/posts\/786\/revisions"}],"predecessor-version":[{"id":1696,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/posts\/786\/revisions\/1696"}],"wp:attachment":[{"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/media?parent=786"}],"wp:term":[{"taxonomy":"category","embeddable":true,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/categories?post=786"},{"taxonomy":"post_tag","embeddable":true,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/tags?post=786"}],"curies":[{"name":"wp","href":"https:\/\/api.w.org\/{rel}","templated":true}]}}