Changeset 1904
- Timestamp:
- Oct 13, 2005, 11:38:46 AM (18 years ago)
- Location:
- inundation/pyvolution
- Files:
-
- 2 edited
Legend:
- Unmodified
- Added
- Removed
-
inundation/pyvolution/least_squares.py
r1903 r1904 920 920 Mandatory input 921 921 time: px1 array of monotonously increasing times (Float) 922 quantities: Dictionary of pxm arrays or 1 pxm array (Float) 922 quantities: Dictionary of arrays or 1 array (Float) 923 The arrays must either have dimensions pxm or mx1. 924 The resulting function will be time dependent in 925 the former case while it will be constan with 926 respect to time in the latter case. 923 927 924 928 Optional input: … … 926 930 vertex_coordinates: mx2 array of coordinates (Float) 927 931 triangles: nx3 array of indices into vertex_coordinates (Int) 928 interpolation_points: array of coordinates to be interpolated to932 interpolation_points: Nx2 array of coordinates to be interpolated to 929 933 verbose: Level of reporting 930 934 … … 1002 1006 self.vertex_coordinates = vertex_coordinates 1003 1007 self.interpolation_points = interpolation_points 1004 self.T = time[:] # Time assumed to be relative to starttime1005 self.index = 0 # Initial time index1008 self.T = time[:] # Time assumed to be relative to starttime 1009 self.index = 0 # Initial time index 1006 1010 self.precomputed_values = {} 1007 1011 … … 1013 1017 raise 'Triangles and vertex_coordinates must be specified' 1014 1018 1015 1016 1019 try: 1017 self.interpolation_points =\ 1018 ensure_numeric(self.interpolation_points) 1020 self.interpolation_points = ensure_numeric(interpolation_points) 1019 1021 except: 1020 1022 msg = 'Interpolation points must be an N x 2 Numeric array '+\ 1021 1023 'or a list of points\n' 1022 msg += 'I got: %s.' %( str(self.interpolation_points)[:60] + '...') 1024 msg += 'I got: %s.' %(str(self.interpolation_points)[:60] +\ 1025 '...') 1023 1026 raise msg 1024 1027 1025 1028 1029 m = len(self.interpolation_points) 1030 p = len(self.T) 1031 1026 1032 for name in quantity_names: 1027 self.precomputed_values[name] =\ 1028 zeros((len(self.T), 1029 len(self.interpolation_points)), 1030 Float) 1033 self.precomputed_values[name] = zeros((p, m), Float) 1031 1034 1032 1035 #Build interpolator 1033 1036 interpol = Interpolation(vertex_coordinates, 1034 1037 triangles, 1035 point_coordinates = self.interpolation_points, 1038 point_coordinates = \ 1039 self.interpolation_points, 1036 1040 alpha = 0, 1037 1041 precrop = False, … … 1039 1043 1040 1044 if verbose: print 'Interpolate' 1041 n = len(self.T)1042 1045 for i, t in enumerate(self.T): 1043 1046 #Interpolate quantities at this timestep 1044 if verbose and i%(( n+10)/10)==0:1045 print ' time step %d of %d' %(i, n)1047 if verbose and i%((p+10)/10)==0: 1048 print ' time step %d of %d' %(i, p) 1046 1049 1047 1050 for name in quantity_names: 1048 self.precomputed_values[name][i, :] =\ 1049 interpol.interpolate(quantities[name][i,:]) 1051 if len(quantities[name].shape) == 2: 1052 result = interpol.interpolate(quantities[name][i,:]) 1053 else: 1054 #Assume no time dependency 1055 result = interpol.interpolate(quantities[name][:]) 1056 1057 self.precomputed_values[name][i, :] = result 1058 1059 1050 1060 1051 1061 #Report -
inundation/pyvolution/test_least_squares.py
r1903 r1904 1129 1129 1130 1130 1131 def test_interpolation_function_spatial_only(self): 1132 """Test spatio-temporal interpolation with constant time 1133 """ 1134 1135 #Three timesteps 1136 time = [1.0, 5.0, 6.0] 1137 1138 1139 #Setup mesh used to represent fitted function 1140 a = [0.0, 0.0] 1141 b = [0.0, 2.0] 1142 c = [2.0, 0.0] 1143 d = [0.0, 4.0] 1144 e = [2.0, 2.0] 1145 f = [4.0, 0.0] 1146 1147 points = [a, b, c, d, e, f] 1148 #bac, bce, ecf, dbe 1149 triangles = [[1,0,2], [1,2,4], [4,2,5], [3,1,4]] 1150 1151 1152 #New datapoints where interpolated values are sought 1153 interpolation_points = [[ 0.0, 0.0], 1154 [ 0.5, 0.5], 1155 [ 0.7, 0.7], 1156 [ 1.0, 0.5], 1157 [ 2.0, 0.4], 1158 [ 2.8, 1.2]] 1159 1160 1161 #One quantity linear in space 1162 Q = linear_function(points) 1163 1164 1165 #Check interpolation of one quantity using interpolaton points 1166 I = Interpolation_function(time, Q, 1167 vertex_coordinates = points, 1168 triangles = triangles, 1169 interpolation_points = interpolation_points, 1170 verbose = False) 1171 1172 1173 answer = linear_function(interpolation_points) 1174 1175 t = time[0] 1176 for j in range(50): #t in [1, 6] 1177 for id in range(len(interpolation_points)): 1178 assert allclose(I(t, id), answer[id]) 1179 1180 t += 0.1 1181 1182 1183 try: 1184 I(1) 1185 except: 1186 pass 1187 else: 1188 raise 'Should raise exception' 1189 1190 1191 1131 1192 def test_interpolation_function(self): 1132 1193 """Test spatio-temporal interpolation … … 1218 1279 1219 1280 1281 1220 1282 def test_fit_and_interpolation_with_different_origins(self): 1221 1283 """Fit a surface to one set of points. Then interpolate that surface
Note: See TracChangeset
for help on using the changeset viewer.