ARTICLE AD BOX
I am working on a 1D rod heat diffusion simulation with a moving heat generation term.
Below is my code for the moving source:
# Functions def wf_pos(t): """ Position of the heat source position in m """ position = np.sqrt(t * dt) # interpolate position_interp = interp1d(t, position) return position_interpI am then integrating the following ODE:
def ode(t, T, Nx, x, dx): """ The ODE from of the spatially discretized PDE """ # laplacian NX = len(T) Tp = np.zeros((NX)) for i in range(1, NX-1): Tp[i]=((T[i-1] - 2 * T[i] + T[i+1])) # finite difference num # compute the heat source at this time # (only need time) x = wf_pos(t) # position of the source (interpolated) lindex = lindex_func(x) # tells us which configuration number we are at Q_val = Qi_func(lindex) # heat generated by source at that point # Laplacian estimation using finite difference method coef = k / (dx**2) Tp = coef * Tp # source term Tp += Q_val # BCs Tp[0] = Tp[1] Tp[NX-1] = Tp[NX-2] # divide by rho*Cp Tp /= rhoCp return Tp # that's what we'll integrateusing
args = (Nx, x, dx) # time integration sol = integrate.solve_ivp( ode, (t_start, t_end), T0, args=args, t_eval=t_eval, rtol=1e-6, atol=1e-6, method='RK45', )My issue is the following. I want a fixed spatial step size Nx = 100, but when plotting the temperature along the rod at different snapshots (5 different snapshots of 3 minutes each):
# plotting for k in range(1, 5): cmap = "RdBu_r" # handle indexation wf_start = wf_pos(t[int((k-1)*t_seg)]/dt) # where the source starts during this time segment wf_end = wf_pos(t[int((k)*t_seg)]/dt) # where the source ends during this time segment print(t[int((k-1)*t_seg) : int((k)*t_seg)].size) print(x_points.size) # interpolate # x_points = interp1d(t[int((k-1)*t_seg) : int((k)*t_seg)] /dt, x_points) print(wf_start) # position of start of wf print(wf_end) # position of end of wf print(np.where(x_points == wf_start)) print(np.where(x_points == wf_end)) # print(np.where(x_points == wf_start)[0].dtype) # print(np.where(x_points == wf_end)[0].dtype) idx_start = int(np.where(x_points == wf_start)[0]) idx_end = int(np.where(x_points == wf_end)[0]) print("idx_start", idx_start) print("idx_end", idx_end) print(x_points[idx_start : idx_end]) print(u[int((k-1)*t_seg) : int((k)*t_seg)]) # plotting # ax[k].plot(x_points[idx_start : idx_end], u[int(k*t_seg)/dt : int((k+1)*t_seg)/dt], label=r'$T_{\mathrm{bed}}(x)$', color='C3', linewidth=1.5) ax[k-1].plot(x_points[idx_start : idx_end], u[idx_start : idx_end], label=r'$T_{\mathrm{bed}}(x)$', color='C3', linewidth=1.5) ax[k-1].set_ylabel(r'$T_{\mathrm{bed}}$ (°C)') ax[k-1].set_xlabel(r'Position along the bed (x)') ax[k-1].set_title('Bed temperature as a function of position for ' + str(int((k-1)*t_seg)) + " to " + str(int((k)*t_seg)) + " min") ax[k-1].legend() ax[k-1].grid(True) fig.tight_layout() plt.show()I'm noticing that because of the position function wf_pos(t) , the position returned at the start and end of different snapshots are not well defined with respect to the x_point spatial vector, defined as
Nx = 100 Lx = 0.15 # m dx = Lx / Nx x_points = np.arange(0, Lx, dx)since it's a square root. I want to be able to handle any function for position as a function of time, which is why I'm trying to use interpolation, but my approach does not work because I'm interpolating two arrays of different shapes (i.e.: x_points has dimension 100, but t has dimension 9000).
It's hard to troubleshoot the error because I am getting TypeError: object of type 'float' has no len(), which is referring to the integration process so I'm not sure when it's happening exactly.
