4

Yes, but you have to mean symplectic on a higher-dimensional phase space than your original problem that includes previous steps too. As I understand there are also some subtle stability issues too. Rather than try to summarize, I'll just refer you to chapter 15, section 4 of Hairer, Lubich, Wanner, Geometric Numerical Integration. That book is a must-have ...


4

You can convert this equation into a Poisson equation and use standard numerical methods to solve it like finite difference to obtain $\theta(x,y)$. If you take a divergence from both side: $$\nabla^{2} \theta = \frac{1}{k}\nabla \cdot \vec{v}$$ Since, you know the $\vec{v}$ you can approximate its divergence as: $$\nabla \cdot \vec{v} = \frac{\partial v_{...


3

The local error in each integration step is composed of 3 parts: the discretization error of the method the floating point error of actually performing the method in precision-limited data types, and specifically for implicit methods, the error in solving the implicit equation. The first error is of size $O(h^{p+1})$, the second of size $\sim \mu$, $\mu$ ...


3

The parameter can be any type, so here I pass in a time-dependent function for p and use it in the differential equation: # Packages library(tidyverse) library(diffeqr) library(JuliaCall) diffeq_setup() # Drift function f <- function(u,p,t){ du1 = p(t) return(c(du1)) } # Diffusion function g <- function(u,p,t){ du1 = 0 # note that there is ...


3

The situation seems to be: You have some input function $x$ which more or less follows a model $\dot x=-ax+bu$. There may be noise involved, so the values of $x$ are not exact, and simply computing difference quotients will in general not be close to the right side of the differential equation. To filter out the noise some averaging is required. This you do ...


2

This doesn't specifically answer your question other than to agree that a finite difference method is an acceptable approach to this problem. Instead, I'm going to propose an alternate approach assuming your main objective is to get a solution to your equation rather than experiment with numerical methods. Particularly for challenging problems like this one,...


2

I can't speak to the use of the FD approach here, but I don't see any reason it wouldn't work. Its fine to have those square terms and ugliness, it just indicates that your governing equations are nonlinear. One way to approach this to say you have a residual, R, that is expressed like so: $$R(u) = a\left(\frac{du}{dz}\right)^2 + b\frac{du}{dz}\frac{d^2u}{...


2

The usual trick is to add more variables that represent the successive derivatives, as in the equations of motion pf physics written as a set of first order ODE of "2" variables: $$\begin{align}\dot{x}&=v\\ \dot{v}&=a(x,t)\end{align}$$ instead of a second order ODE of 1 variable $$\ddot{x}=a(x,t) \ .$$ The dot represents the time derivative. So, ...


2

In fact, your equation is a non-linear advection-diffusion. Due to the fact that your problem is time-dependent, it could be easily solved by finite-difference: $$\frac{\partial z}{\partial t} = -C (\sin(\omega t))^{m} x^{hm} n (\frac{\partial z}{\partial x})^{n-1} \frac{\partial^{2} z}{\partial x^{2}} -C (\sin(\omega t))^{m} x^{hm} h m x^{-1} (\frac{\...


2

You will certainly benefit from parallelization of your code (or using/switching libraries that offer a native parallelization). So it might be worth looking into parallelization depending on your priorities in terms of time and possibility of new upcoming simulations. In terms of external computing resources, I would say, that you are looking for a very ...


1

I think the problem with initial conditions is more related with the sum of very small numbers with other that are only small (I suspect the problem resides in terms like gam^8, since gam is small). The resulting numeric inaccuracy will bring problems. I've played around with other initial conditions and integration times and in octave (that is similar but ...


1

Take the first expression and start to reduce the $x$ values to $x_i$, \begin{align} \frac{u_{i+1}^n - u_i^n}{x_{i+\frac{1}{2}}} - \frac{u_i^n - u_{i-1}^n}{x_{i-\frac{1}{2}}} &= \frac{(x_i-\frac12Δx)(u_{i+1}^n - u_i^n) - (x_i+\frac12Δx)(u_i^n - u_{i-1}^n)}{x_{i-\frac{1}{2}}x_{i+\frac{1}{2}}} \\ &=\frac{x_i}{x_i^2-\frac14Δx^2}(u_{i+1}^n - 2u_i^n + u_{...


1

If I understand your problem correctly, you just need to do parameter identification from a known model and non-noisy data. The standard way to do this is via a nonlinear least-squares framework. To do this, after solving the ODE for your given choice of $\rho$ and $f$ at some number of time points, you create a cost function that is a function of $\rho$ and ...


1

First of all, you need to rewrite your non-linear differential equation in a set of first order differential equations since the solve_ivp routine solves problems of the form $u'(x) = f(x,u), \ u(0) = u_{0}$ where $u(x)$ can be a vector. So define $u_{1}(x) = y(x)$ and $u_{2}(x) = y'(x)$ and you can write your original equation as \begin{array}{lll} u_{1}'(x)...


1

This comes from von Neumann stability analysis. You discretize the partial differential equation and look at the error equation, which is the difference between the exact solution to the finite difference approximation and the actual equation. You need the error to shrink in time, and if you assume that the error behaves like a fourier series, you can ...


1

I'm going to turn my comment into an answer. The error order tells you up to which order (exclusive) the discrete solution corresponds to the exact solution. An $\mathcal{O}(\Delta t^4)$ method reproduces the orders 0-3 of the Taylor expansion of the function around $\Delta t$ (i.e. in the expression $f(t+\Delta t) = f(t) +\dots$). Put bluntly (and ...


1

For the interpretation where the left and right "virtual" compartments are infinite sinks at level Z0=Z9=0, and the walls between the compartments j and j+1 are removed at time ts[j] where also the external filling of compartment j starts, I get the shortened implementation def ODE_func(Z,t,Q,Zd,Z0,dts,ts,cf,t_copy): # Declare some boundary conditions: ...


Only top voted, non community-wiki answers of a minimum length are eligible