Approximate Fourier series recursion for problems involving temporal fractional calculus

Fractional calculus constitutive models are of wide applicability. They have been used, for example, for: xanthan gum, bread dough, synthetic polymers (e.g. nylon), single collagen fibrils, semi-hard zero-fat cheese, pure ice at -35 degrees C, asphalt sealants, epithelial cells, gels, polymers, concrete, asphalt, rock mass, waxy crude oil, breast tissue cells, lung parenchyma, and red blood cell membranes (e.g. Bonfanti et al. (2020)). In this note we are motivated by applications to dispersive media such as lossy dielectrics and viscoelastic polymers, both of which involve constitutive laws with fading memory convolution Volterra integrals, and both of which have been modelled by fractional calculus - which implies a weakly singular Volterra kernel. To step forward in time with such a model the entire hereditary Volterra integral needs to be numerically re-computed at each time step. Done naively with standard quadrature this involves O(N2) operations and O(N) storage to compute over N time steps, and this severely limits the usefulness of such simple methods in large scale computations. Several alternatives addressing these shortcomings of this naive approach have been offered in the literature, all of which can provide remedies. Here we propose another method, but ours is rather different in that we use a Fourier series proxy over much of the integration range. This allows for a recursive update to the 'history integral', with a much smaller effort expended in standard quadrature near the singularity. The recursive update is basically the same as that used when the Volterra kernel comprises decaying exponentials, a method which is already employed in some commercial codes. Given that this partial implementation is already in use, we intend that our method therefore has a lower 'barrier to entry' for incorporating fractional calculus models into commercial codes: it requires only minor surgery on the weakly singular kernel in order to apply the Fourier series proxy. This study is presented in the spirit of a discussion piece - there are several directions one could take from the basic observations presented here. We have aimed for simplicity and economy of presentation, and included just enough numerical results to give evidence that the method works. In summary, our method mitigates the O(N) storage and O(N2) storage issues using an easily implemented proxy. It can also be implemented for variable time steps. The code is available from https://github.com/variatio nalform/fouvol and https://hub.docker.com/u/variationalform.