Equation is not affine w.r.t the target,

[gfkdliucheng]

def solver(I, w, dt, T): """ Solve u'' + w2u = 0 for t in (0,T], u(0)=I and u'(0)=0, by a central finite difference method with time step dt. """ dt = float(dt) Nt = int(round(T/dt)) t = Dimension('t', spacing=Constant('h_t')) u = TimeFunction(name='u', dimensions=(t,), shape=(Nt+1,), space_order=2) u.data[:] = I eqn = u.dt2 + (w2)u stencil = Eq(u.forward, solve(eqn, u.forward)) op = Operator(stencil) op.apply(h_t=dt, t_M=Nt-1) return u.data, np.linspace(0, Nt*dt, Nt+1)

I = 1 w = 2np.pi dt = 0.05 num_periods = 5 P = 2np.pi/w # one period T = P*num_periods u, t = solver(I, w, dt, T)

errors are as follows： Equation is not affine w.r.t the target, falling back to standardsympy.solve that may be slow /tmp/devito-jitcache-uid1000/d8ebcf7bd2f4100cee32a40cc01b7a94b4c2516d.c: In function 'Kernel': /tmp/devito-jitcache-uid1000/d8ebcf7bd2f4100cee32a40cc01b7a94b4c2516d.c:41:57: error: array subscript is not an integer 41 | u[t + 1] = -3.94784176043574e+1Fr0h_th_t + 2.0Fu[-39r0] - 1.0Fu[-39r0 - 1]; | ^ /tmp/devito-jitcache-uid1000/d8ebcf7bd2f4100cee32a40cc01b7a94b4c2516d.c:41:74: error: array subscript is not an integer 41 | u[t + 1] = -3.94784176043574e+1Fr0h_th_t + 2.0Fu[-39r0] - 1.0Fu[-39r0 - 1]; | ^ FAILED compiler invocation: gcc -march=native -O3 -g -fPIC -Wall -std=c99 -Wno-unused-result -Wno-unused-variable -Wno-unused-but-set-variable -ffast-math -mprefer-vector-width=512 -shared -fopenmp /tmp/devito-jitcache-uid1000/d8ebcf7bd2f4100cee32a40cc01b7a94b4c2516d.c -lm -o /tmp/devito-jitcache-uid1000/d8ebcf7bd2f4100cee32a40cc01b7a94b4c2516d.so

[mloubout]

t = Dimension('t', spacing=Constant('h_t'))
u = TimeFunction(name='u', dimensions=(t,),
                 shape=(Nt+1,), space_order=2)

• Since you are defining a TimeFunction with a single dimension the dimension needs to be a TImeDimension for the function to be valid or it will behave weirdly:
• You need at least time_order=2 for .dt2 to be properly defined (and with the proper time dimensions above)

    t = TimeDimension('t', spacing=Constant('h_t'))

    u = TimeFunction(name='u', dimensions=(t,),
                    shape=(Nt+1,), space_order=2, time_order=2)

and should be all good

[gfkdliucheng]

e it's not a program issue, this is from the official dem

The problem still exists, it seems like it's not a program issue, this is the official demo.

[mloubout]

The problem still exists

Not sure what you mean, if you follow the above change it all works fine

In [16]: from devito import *
    ...: import numpy as np
    ...: 
    ...: def solver(I, w, dt, T):
    ...:     """
    ...:     Solve u'' + w**2*u = 0 for t in (0,T], u(0)=I and u'(0)=0,
    ...:     by a central finite difference method with time step dt.
    ...:     """
    ...:     dt = float(dt)
    ...:     Nt = int(round(T/dt))
    ...:     t = TimeDimension('t', spacing=Constant('h_t'))
    ...: 
    ...:     u = TimeFunction(name='u', dimensions=(t,),
    ...:                     shape=(Nt+1,), space_order=2, time_order=2)
    ...: 
    ...:     u.data[:] = I
    ...:     eqn = u.dt2 + (w**2)*u
    ...:     stencil = Eq(u.forward, solve(eqn, u.forward))
    ...:     op = Operator(stencil)
    ...:     op.apply(h_t=dt, t_M=Nt-1)
    ...:     return u.data, np.linspace(0, Nt*dt, Nt+1)
    ...: 
    ...: I = 1
    ...: w = 2*np.pi
    ...: dt = 0.05
    ...: num_periods = 5
    ...: P = 2*np.pi/w # one period
    ...: T = P*num_periods
    ...: u, t = solver(I, w, dt, T)
Operator `Kernel` ran in 0.01 s

[gfkdliucheng]

The problem still exists

Not sure what you mean, if you follow the above change it all works fine

In [16]: from devito import *
    ...: import numpy as np
    ...: 
    ...: def solver(I, w, dt, T):
    ...:     """
    ...:     Solve u'' + w**2*u = 0 for t in (0,T], u(0)=I and u'(0)=0,
    ...:     by a central finite difference method with time step dt.
    ...:     """
    ...:     dt = float(dt)
    ...:     Nt = int(round(T/dt))
    ...:     t = TimeDimension('t', spacing=Constant('h_t'))
    ...: 
    ...:     u = TimeFunction(name='u', dimensions=(t,),
    ...:                     shape=(Nt+1,), space_order=2, time_order=2)
    ...: 
    ...:     u.data[:] = I
    ...:     eqn = u.dt2 + (w**2)*u
    ...:     stencil = Eq(u.forward, solve(eqn, u.forward))
    ...:     op = Operator(stencil)
    ...:     op.apply(h_t=dt, t_M=Nt-1)
    ...:     return u.data, np.linspace(0, Nt*dt, Nt+1)
    ...: 
    ...: I = 1
    ...: w = 2*np.pi
    ...: dt = 0.05
    ...: num_periods = 5
    ...: P = 2*np.pi/w # one period
    ...: T = P*num_periods
    ...: u, t = solver(I, w, dt, T)
Operator `Kernel` ran in 0.01 s

Thank you, bro