SHARE
TWEET

Untitled

a guest Feb 17th, 2020 73 Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
  1. import math as m
  2. import matplotlib.pyplot as plt
  3.  
  4. u0 = 0.1
  5. C = 30
  6.  
  7. def f1(s):
  8.     return s**2
  9.  
  10. def f2(s):
  11.     return (1-s)**2
  12.  
  13. def fi(s, c):
  14.     a = mu(c) * ((s - 0.2) / 0.81)**3
  15.     b = a + ((0.74 - s) / 0.715)**3
  16.     return  a / b
  17.  
  18. def mu(c):
  19.     return 20 * (1 + 5*c)
  20.  
  21. def a(c):
  22.     return 0.1*c
  23.  
  24. tau = 0.0002
  25. h = 0.01
  26.  
  27. N = 100
  28. T = 1000
  29. m = 0.2
  30. U = 1
  31. D = 0.0001
  32.  
  33. c = [[None for i in range(N)] for j in range(T)]
  34. s = [[None for i in range(N)] for j in range(T)]
  35.  
  36. for i in range (N):
  37.     s[0][i] = 0.74
  38.     c[0][i] = 1
  39.  
  40. for j in range (T-1):
  41.     s[j][0] = 0.2
  42.     c[j][0] = 0
  43.     c[j][N-1] = 1
  44.     for i in range(1, N-1):
  45.         s[j+1][i] = s[j][i] - (fi(s[j][i], c[j][i]) - fi(s[j][i-1], c[j][i-1]))/h*tau/m
  46.         c1 = s[j][i]*c[j][i] + 0.1*c[j][i]/m + (-U*(fi(s[j][i], c[j][i]) - fi(s[j][i-1], c[j][i])*c[j][i-1]) + D*(c[j][i-1] - 2*c[j][i] + c[j][i+1])/h**2)*tau/m
  47.         c[j+1][i] = c1/(s[j+1][i] + 0.1/m)
  48.  
  49. print(s[300])
  50. print(c[300])
RAW Paste Data
We use cookies for various purposes including analytics. By continuing to use Pastebin, you agree to our use of cookies as described in the Cookies Policy. OK, I Understand
Top