Skip to content

Commit

Permalink
test
Browse files Browse the repository at this point in the history
  • Loading branch information
statespacedev committed Sep 21, 2024
1 parent 9f566f2 commit d08db96
Show file tree
Hide file tree
Showing 3 changed files with 54 additions and 57 deletions.
18 changes: 9 additions & 9 deletions pdp10/algebra/ipwr.for
Original file line number Diff line number Diff line change
@@ -1,10 +1,10 @@

program hello
integer :: i
real :: pi
pi = 3.1415927
i = 0
i = i + 1
print *, 'Hello, World!', pi
end program hello

program hello
integer i
real pi
pi = 3.1415927
i = 0
i = i + 1
print *, 'Hello, World!', pi
end program hello

40 changes: 20 additions & 20 deletions pdp10/algebra/ipwr.py
Original file line number Diff line number Diff line change
@@ -1,21 +1,21 @@
""" imaginary powers of base 10. starting with 10**(i/1024), and squaring successively ten times. this matches with feynman's table 22-3"""
from math import sqrt

# non rounded version
y = .00225 # should be .0022486
x = sqrt(1. - y**2)
for _ in range(11):
print('%10.7f %10.7f' % (x, y))
x2 = x**2 - y**2
y2 = 2 * x * y
x, y = x2, y2

# rounded version that more closely matches feynman's table 22-3
y = .00225 # should be .0022486
x = round(sqrt(1. - y**2), 7)
for _ in range(11):
print('%10.5f %10.5f' % (x, y))
x2 = round(x**2 - y**2, 7)
y2 = round(2 * x * y, 7)
x, y = x2, y2
""" imaginary powers of base 10. starting with 10**(i/1024), and squaring successively ten times. this matches with feynman's table 22-3"""
from math import sqrt

# non rounded version
y = .00225 # should be .0022486
x = sqrt(1. - y**2)
for _ in range(11):
print('%10.7f %10.7f' % (x, y))
x2 = x**2 - y**2
y2 = 2 * x * y
x, y = x2, y2

# rounded version that more closely matches feynman's table 22-3
y = .00225 # should be .0022486
x = round(sqrt(1. - y**2), 7)
for _ in range(11):
print('%10.5f %10.5f' % (x, y))
x2 = round(x**2 - y**2, 7)
y2 = round(2 * x * y, 7)
x, y = x2, y2

53 changes: 25 additions & 28 deletions pdp10/algebra/simplex.py
Original file line number Diff line number Diff line change
@@ -1,28 +1,25 @@
"""simplest possible simplex implementation"""
import numpy as np


def main():
A = np.array([[-2, 1, 1, 0, 0], [-1, 2, 0, 1, 0], [1, 0, 0, 0, 1]])
b = np.array([2, 7, 3])
c = np.array([-1, -2, 0, 0, 0])
x = simplex(A, b, c)
return x


def simplex(A, b, c):
"""revised simplex method. """
m, n, o = A.shape[0], A.shape[1], A.shape[1] - A.shape[0]
ns, v1, v2 = [i for i in range(n)], np.array(c[o:]), np.array(c[:o])
Binv = np.linalg.inv(A[:, ns[o:]])
while not np.min(v2 - v1 @ Binv @ A[:, ns[:o]]) > 0:
n1 = np.argmin(v2 - v1 @ Binv @ A[:, ns[:o]])
t1, t2 = Binv @ b, Binv @ A[:, ns[n1]]
n2 = np.argmin([t1[i] / t2[i] if t2[i] > 0 else np.inf for i in range(m)])
ns[n1], ns[n2 + o], v1[n2], v2[n1] = ns[n2 + o], ns[n1], v2[n1], v1[n2]
Binv = np.linalg.inv(A[:, ns[o:]])
return Binv @ b


if __name__ == "__main__":
main()
"""simplest possible simplex implementation"""
import numpy as np

def main():
A = np.array([[-2, 1, 1, 0, 0], [-1, 2, 0, 1, 0], [1, 0, 0, 0, 1]])
b = np.array([2, 7, 3])
c = np.array([-1, -2, 0, 0, 0])
x = simplex(A, b, c)
return x

def simplex(A, b, c):
"""revised simplex method. """
m, n, o = A.shape[0], A.shape[1], A.shape[1] - A.shape[0]
ns, v1, v2 = [i for i in range(n)], np.array(c[o:]), np.array(c[:o])
Binv = np.linalg.inv(A[:, ns[o:]])
while not np.min(v2 - v1 @ Binv @ A[:, ns[:o]]) > 0:
n1 = np.argmin(v2 - v1 @ Binv @ A[:, ns[:o]])
t1, t2 = Binv @ b, Binv @ A[:, ns[n1]]
n2 = np.argmin([t1[i] / t2[i] if t2[i] > 0 else np.inf for i in range(m)])
ns[n1], ns[n2 + o], v1[n2], v2[n1] = ns[n2 + o], ns[n1], v2[n1], v1[n2]
Binv = np.linalg.inv(A[:, ns[o:]])
return Binv @ b

if __name__ == "__main__":
main()

0 comments on commit d08db96

Please sign in to comment.