 Timestamp:
 Jun 22, 2010, 5:30:32 PM (13 years ago)
 File:

 1 edited
Legend:
 Unmodified
 Added
 Removed

anuga_work/development/2010projects/anuga_1d/base/quantity.py
r7855 r7868 70 70 #flux calculations and forcing functions 71 71 72 N = domain.number_of_elements 72 self.N = domain.number_of_elements 73 N = self.N 73 74 self.explicit_update = numpy.zeros(N, numpy.float ) 74 75 self.semi_implicit_update = numpy.zeros(N, numpy.float ) … … 78 79 self.qmin = numpy.zeros(self.centroid_values.shape, numpy.float) 79 80 80 self.beta = domain.beta 81 #These are taken from domain but can be set for each quantity 82 # if so desired 83 self.beta = 0.0 84 self.limiter = 'vanleer' 81 85 82 86 … … 604 608 return integral 605 609 610 def get_beta(self,beta): 611 """Get limiting parameter 612 """ 613 614 return self.beta 615 616 def set_beta(self,beta): 617 """Set limiting parameter 618 """ 619 620 #Probably should test that it is not too large 621 self.beta = beta 622 623 624 def get_limiter(self): 625 return self.limiter 626 627 def set_limiter(self,limiter): 628 629 possible_limiters = \ 630 ['pyvolution', 'minmod_steve', 'minmod', 'minmod_kurganov', 'superbee', 'vanleer', 'vanalbada'] 631 632 if limiter in possible_limiters: 633 self.limiter = limiter 634 return 635 636 msg = '%s is an incorrect limiter type.\n'% limiter 637 msg += 'Possible types are: '+ ", ".join(["%s" % el for el in possible_limiters]) 638 raise Exception, msg 606 639 607 640 def update(self, timestep): … … 795 828 796 829 #Gradient 797 G[1:1] = xmic( self. domain.beta, d1, d2 )830 G[1:1] = xmic( self.beta, d1, d2 ) 798 831 799 832 … … 818 851 """ 819 852 820 if self. domain.limiter == "pyvolution":853 if self.limiter == "pyvolution": 821 854 self.limit_pyvolution() 822 855 823 elif self. domain.limiter == "minmod_steve":856 elif self.limiter == "minmod_steve": 824 857 self.limit_minmod() 825 858 … … 985 1018 from limiters_python import minmod, minmod_kurganov, minmod_kurganov_old, maxmod, vanleer, vanalbada 986 1019 987 limiter = self. domain.limiter1020 limiter = self.get_limiter() 988 1021 #print limiter 989 1022 990 1023 #print 'limit_range' 991 N = self. domain.number_of_elements1024 N = self.N 992 1025 qc = self.centroid_values 993 1026 qv = self.vertex_values … … 1013 1046 phi[N1] = (qc[N1]  qc[N2])/(xc[N1]  xc[N2]) 1014 1047 1048 1015 1049 if limiter == "minmod": 1016 1050 phi[1:1] = minmod(beta_p[1:1],beta_m[1:1]) … … 1023 1057 1024 1058 elif limiter == "minmod_kurganov": 1025 theta = 2.01059 theta = self.beta 1026 1060 phi[1:1] = minmod_kurganov(theta*beta_p[1:1],theta*beta_m[1:1], beta_x[1:1]) 1027 1061 … … 1034 1068 msg = 'Unknown limiter' 1035 1069 raise Exception, msg 1036 1070 1071 1072 1037 1073 qv[:,0] = qc + phi*dx[:,0] 1038 1074 qv[:,1] = qc + phi*dx[:,1] 1039 1075 1040 1076 1041 #Phi limiter 1042 # for k in range(N): 1043 # n0 = self.domain.neighbours[k,0] 1044 # n1 = self.domain.neighbours[k,1] 1045 # if n0 < 0: 1046 # phi = (qc[k+1]  qc[k])/(xc[k+1]  xc[k]) 1047 # elif n1 < 0: 1048 # phi = (qc[k]  qc[k1])/(xc[k]  xc[k1]) 1049 # #elif (self.domain.wet_nodes[k,0] == 2) & (self.domain.wet_nodes[k,1] == 2): 1050 # # phi = 0.0 1051 # else: 1052 # if limiter == "minmod": 1053 # phi = minmod(beta_p[k],beta_m[k]) 1054 # 1055 # elif limiter == "minmod_kurganov":#Change this 1056 # # Also known as monotonized central difference limiter 1057 # # if theta = 2.0 1058 # theta = 2.0 1059 # phi = minmod_kurganov(theta*beta_p[k],theta*beta_m[k],beta_x[k]) 1060 # 1061 # elif limiter == "superbee": 1062 # slope1 = minmod(beta_m[k],2.0*beta_p[k]) 1063 # slope2 = minmod(2.0*beta_m[k],beta_p[k]) 1064 # phi = maxmod(slope1,slope2) 1065 # 1066 # elif limiter == "vanleer": 1067 # phi = vanleer(beta_p[k],beta_m[k]) 1068 # 1069 # elif limiter == "vanalbada": 1070 # phi = vanalbada(beta_m[k],beta_p[k]) 1071 # 1072 # for i in range(2): 1073 # qv[k,i] = qc[k] + phi*dx[k,i] 1077 1074 1078 1075 1079 def limit_steve_slope(self):
Note: See TracChangeset
for help on using the changeset viewer.