I wrote this function called cfd() as a central finite-difference algorithm.
The input is a 1-D array ‘x’.
The for-loop starts with creating x_fwd and x_bwd as copies of the original input x. My problem is right after that when I reassign the i-th element of x_fwd to a new value (x_fwd[i] = x_fwd[i] + (h/2)). After this line of code, the print statements reveal that the i-th element of both x and x_bwd have also been given a new value, even though the only code executed was for x_fwd. Why??
See code below.
I need x to stay the same always and I need to change elements of x_bwd and x_fwd independently of each other.
# central finite difference function
def cfd(x):
print("entering cfd()")
jacobian = np.zeros_like(x)
print("initial jac: ", jacobian)
h = 0.2
print("h>> ", h)
for i in range(0, 10):
print("\niteration: ", i)
print("x>> ", x)
x_fwd = x # x-forward as a copy of the input x
x_bwd = x # x-backward as a copy of the input x
print("x_fwd>> ", x_fwd)
print("x_bwd>> ", x_bwd)
print("Add h/2 to x_fwd[i]:")
x_fwd[i] = x_fwd[i] + (h/2) # <--- add h/2 to the i-th element of x_fwd
print("x_fwd>> ", x_fwd) # <--- this and the next two print lines reveal that the
print("x_bwd>> ", x_bwd) # i-th element of x_bwd and x are also changed. Why??
print("x>> ", x)
print("Subtract h/2 from x_bwd[i]:")
x_bwd[i] = x_bwd[i] - (h/2)
print("x_fwd>> ", x_fwd)
print("x_bwd>> ", x_bwd)
print("x>> ", x)
print("f at x_fwd and x_bwd:")
print("f(x_fwd)>> ", f(x_fwd))
print("f(x_bwd)>> ", f(x_bwd))
# the i-th element of the jacobian calcualted by central finite-difference
jacobian[i] = (f(x_fwd) - f(x_bwd)) / h
print("new jacobian: ", jacobian)
return jacobian
x0 = [0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2]
cfd(x0)
You believe that is what it is doing, but it is not doing that
x_fwd and x_bwd are references to the same object that x points to; if you are familiar with programming languages which offer pointers, imagine that x is a pointer and you’ve just copied the pointer (twice), but not the thing it points to.
If you want to copy x, you’ll need to use something like this.
x_fwd = x # x-forward as a copy of the input x
x_bwd = x # x-backward as a copy of the input x
These are not copies. They make the names x_fwd and xbwd refer to
the same array that x does. Thus you see the changes via any of
these names, because there’s only one array in use.
You probably want:
x_fwd = x.copy() # x-forward as a copy of the input x
x_bwd = x.copy() # x-backward as a copy of the input x
Just to note, x there is actually a Python list, not a NumPy array, but it also has a .copy() method so your solution still works the same.
Note, though, that if the list contained other non-immutable objects, (e.g. nested lists, dictionaries, arrays, custom class instances, etc), those objects would not be copied recursively and would point to the same object in both copies, unless you used copy.deepcopy(). Not an issue in your case if x is required to contain only numeric types, but something to keep in mind in case you run into this again in another context.
Right, that would make the most logical sense (as jacobian is an array, and that would allow solving the problem much more efficiently and concisely than a for loop via vectorization), but cfd is in fact passed a list for x, and is never cast to an array:
Converting x to an array should hopefully allow solving the problem more cleanly and efficiently using vectorization, and also both inherently make a (shallow) copy of the original list, and avoid mutating it in place to begin with. However, it is not clear what function the cryptically-named f is nor what it should do, so I cannot be sure, nor produce an example to benchmark.
