Last change
on this file since 7837 was
7837,
checked in by mungkasi, 13 years ago

Again, adding some codes for 1d problems on debris avalanche and periodic waves.

File size:
1.1 KB

Line  

1  ## module gaussPivot 

2  ''' x = gaussPivot(a,b,tol=1.0e9). 

3  Solves [a]{x} = {b} by Gauss elimination with 

4  scaled row pivoting 

5  ''' 

6  from Numeric import * 

7  import swap 

8  import error 

9  

10  def gaussPivot(a,b,tol=1.0e12): 

11  n = len(b) 

12  

13  # Set up scale factors 

14  s = zeros((n),Float) 

15  for i in range(n): 

16  s[i] = max(abs(a[i,:])) 

17  

18  for k in range(0,n1): 

19  

20  # Row interchange, if needed 

21  p = int(argmax(abs(a[k:n,k])/s[k:n])) + k 

22  if abs(a[p,k]) < tol: error.err('Matrix is singular') 

23  if p != k: 

24  swap.swapRows(b,k,p) 

25  swap.swapRows(s,k,p) 

26  swap.swapRows(a,k,p) 

27  

28  # Elimination 

29  for i in range(k+1,n): 

30  if a[i,k] != 0.0: 

31  lam = a[i,k]/a[k,k] 

32  a[i,k+1:n] = a [i,k+1:n]  lam*a[k,k+1:n] 

33  b[i] = b[i]  lam*b[k] 

34  if abs(a[n1,n1]) < tol: error.err('Matrix is singular') 

35  

36  # Back substitution 

37  for k in range(n1,1,1): 

38  b[k] = (b[k]  dot(a[k,k+1:n],b[k+1:n]))/a[k,k] 

39  return b 

40  

41  

42  

43  

Note: See
TracBrowser
for help on using the repository browser.