Figure A.6:

Relative error in the effectiveness factor versus number of collocation points.

Code for Figure A.6

Text of the GNU GPL.

main.py

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
# Converted from etaerr.m - effectiveness factor collocation error vs npts
# [makes] etaerr.dat
import numpy as np
from scipy.optimize import fsolve
from misc import save_ascii, colloc

Phi = 10.0; n = 1
etaan = 1.0/Phi * (1.0/np.tanh(3*Phi) - 1.0/(3*Phi))

colvec    = [5, 10, 15, 20, 25, 30]
etaerror  = np.zeros(len(colvec))

for ii, npts in enumerate(colvec):
    R0, A, B, Q = colloc(npts - 2, 'left', 'right')
    R = R0 * 3.0;  A = A/3.0;  B = B/9.0;  Q = Q*3.0

    def pellet(c):
        res = B @ c + 2.0 * A @ c / R - 2.0/(n+1) * Phi**2 * c**n
        res[0]  = A[0, :] @ c
        res[-1] = 1.0 - c[-1]
        return res

    c0 = 0.5 * np.linspace(0, 1, npts)
    c, _, ier, _ = fsolve(pellet, c0, full_output=True)
    eta = (n+1)/2.0 * (A[npts-1, :] @ c) / Phi**2
    etaerror[ii] = abs(eta - etaan) / etaan

table = np.column_stack([colvec, etaerror])
save_ascii('etaerr.dat', table)