Add Z dimension to 2D Homogeneous Transform Matrix

This is something I've been trying to work out for a while - I'm trying to go from a 2D x,y,w homogeneous transform matrix (eg. Android's graphics.matrix) and turn it into a 3D x,y,z,w transform matrix (eg. Android's opengl.matrix)

Convert an android.graphics.Matrix to a GL mat4? seems to be a similar question but the answer only applies to affine transformations while I need to handle perspective transformations in x and y. This also seems like something more general than a Android specific implementation.

I've tried making a system of equations to solve for what the 3D matrix values need to be such that for five points, the results for X and Y equal their value after the 2D transformation is applied and Z always equals one, but if the Z value is always left the same (what I want to happen) then the matrix is not invertible.

import numpy as np

# arbitrary points to sample
# need 5 points to have 5 * 3d = 15 equations to solve for 15 matrix coefficents

x0 = -1
y0 = 1
z0 = 1
u0 = -1
v0 = 1
w0 = 1

x1 = 1
y1 = 1
z1 = 1
u1 = 1
v1 = 1
w1 = 1

x2 = 1
y2 = -1
z2 = 1
u2 = 1
v2 = -1
w2 = 1

x3 = -1
y3 = -1
z3 = 1
u3 = -1
v3 = -1
w3 = 1

x4 = 0
y4 = 0
z4 = 1
u4 = 0
v4 = 0
w4 = 1

# The following numbers seem to allow the calculation to finish but are not what i want
# x4 = 21
# y4 = 15
# z4 = 36
# u4 = 12
# v4 = 21
# w4 = 31

# matrix made by extending https://stackoverflow.com/a/57280136 into 3d

A = np.matrix([
    [x0, y0, z0, 1, 0, 0, 0, 0, 0, 0, 0, 0, -x0*u0, -y0*u0, -z0*u0],
    [x1, y1, z1, 1, 0, 0, 0, 0, 0, 0, 0, 0, -x1*u1, -y1*u1, -z1*u1],
    [x2, y2, z2, 1, 0, 0, 0, 0, 0, 0, 0, 0, -x2*u2, -y2*u2, -z2*u2],
    [x3, y3, z3, 1, 0, 0, 0, 0, 0, 0, 0, 0, -x3*u3, -y3*u3, -z3*u3],
    [x4, y4, z4, 1, 0, 0, 0, 0, 0, 0, 0, 0, -x4*u4, -y4*u4, -z4*u4],
    [0, 0, 0, 0, x0, y0, z0, 1, 0, 0, 0, 0, -x0*v0, -y0*v0, -z0*v0],
    [0, 0, 0, 0, x1, y1, z1, 1, 0, 0, 0, 0, -x1*v1, -y1*v1, -z1*v1],
    [0, 0, 0, 0, x2, y2, z2, 1, 0, 0, 0, 0, -x2*v2, -y2*v2, -z2*v2],
    [0, 0, 0, 0, x3, y3, z3, 1, 0, 0, 0, 0, -x3*v3, -y3*v3, -z3*v3],
    [0, 0, 0, 0, x4, y4, z4, 1, 0, 0, 0, 0, -x4*v4, -y4*v4, -z4*v4],
    [0, 0, 0, 0, 0, 0, 0, 0, x0, y0, z0, 1, -x0*w0, -y0*w0, -z0*w0],
    [0, 0, 0, 0, 0, 0, 0, 0, x1, y1, z1, 1, -x1*w1, -y1*w1, -z1*w1],
    [0, 0, 0, 0, 0, 0, 0, 0, x2, y2, z2, 1, -x2*w2, -y2*w2, -z2*w2],
    [0, 0, 0, 0, 0, 0, 0, 0, x3, y3, z3, 1, -x3*w3, -y3*w3, -z3*w3],
    [0, 0, 0, 0, 0, 0, 0, 0, x4, y4, z4, 1, -x4*w4, -y4*w4, -z4*w4]
])

print(A)
print(np.linalg.det(A)) # zero

b = np.array([u0, u1, u2, u3, u4, v0, v1, v2, v3, v4, w0, w1, w2, w3, w4])

c = np.linalg.solve(A, b) # crashes here

mat3d = np.matrix([
    [c[0], c[1], c[2], c[3]],
    [c[4], c[5], c[6], c[7]],
    [c[8], c[9], c[10], c[11]],
    [c[12], c[13], c[14], 1]
])

print(mat3d)

Is there a way to reliably extend a 2D homogeneous coordinates matrix to a 3D one?

Thanks!