Skip to content

Commit

Permalink
Merge 1913336 into 83708b5
Browse files Browse the repository at this point in the history
  • Loading branch information
@PetarMilevski
PetarMilevski committed Nov 19, 2017
2 parents 83708b5 + 1913336 commit 5859071
Show file tree
Hide file tree
Showing 9 changed files with 397 additions and 229 deletions.
2 changes: 1 addition & 1 deletion anuga/parallel/__init__.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -21,7 +21,7 @@
from parallel_shallow_water import Parallel_domain as Parallel_shallow_water_domain
from parallel_advection import Parallel_domain as Parallel_advection_domain
from parallel_operator_factory import Inlet_operator, Boyd_box_operator, Boyd_pipe_operator
from parallel_operator_factory import Weir_orifice_trapezoid_operator #added by PM 22/10/2013
from parallel_operator_factory import Weir_orifice_trapezoid_operator
else:
from anuga import rectangular_cross as parallel_rectangle
from anuga import Domain as Parallel_shallow_water_domain
Expand Down
12 changes: 12 additions & 0 deletions anuga/parallel/parallel_operator_factory.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -462,16 +462,20 @@ def Boyd_pipe_operator(domain,
def Weir_orifice_trapezoid_operator(domain,
losses,
width,
blockage=0.0,
barrels=1.0,
z1=None,
z2=None,
height=None,
end_points=None,
exchange_lines=None,
enquiry_points=None,
invert_elevations=None,
#culvert_slope=None,
apron=0.1,
manning=0.013,
enquiry_gap=0.0,
smoothing_timescale=0.0,
use_momentum_jet=True,
use_velocity_head=True,
description=None,
Expand All @@ -489,15 +493,19 @@ def Weir_orifice_trapezoid_operator(domain,
losses=losses,
width=width,
height=height,
blockage=blockage,
barrels=barrels,
z1=z1,
z2=z2,
#culvert_slope=culvert_slope,
end_points=end_points,
exchange_lines=exchange_lines,
enquiry_points=enquiry_points,
invert_elevations=invert_elevations,
apron=apron,
manning=manning,
enquiry_gap=enquiry_gap,
smoothing_timescale=smoothing_timescale,
use_momentum_jet=use_momentum_jet,
use_velocity_head=use_velocity_head,
description=description,
Expand Down Expand Up @@ -594,15 +602,19 @@ def Weir_orifice_trapezoid_operator(domain,
losses=losses,
width=width,
height=height,
blockage=blockage,
barrels=barrels,
z1=z1,
z2=z2,
#culvert_slope=culvert_slope,
end_points=end_points,
exchange_lines=exchange_lines,
enquiry_points=enquiry_points,
invert_elevations=invert_elevations,
apron=apron,
manning=manning,
enquiry_gap=enquiry_gap,
smoothing_timescale=smoothing_timescale,
use_momentum_jet=use_momentum_jet,
use_velocity_head=use_velocity_head,
description=description,
Expand Down
10 changes: 6 additions & 4 deletions anuga/parallel/parallel_structure_operator.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -597,10 +597,12 @@ def statistics(self):
message += 'Culvert Blockage: %s\n'% self.blockage
message += 'No. of barrels: %s\n'% self.barrels
else:
message += 'Culvert Height: %s\n'% self.height
message += 'Culvert Width: %s\n'% self.width
message += 'Batter Slope 1: %s\n'% self.z1
message += 'Batter Slope 2: %s\n'% self.z2
message += 'Culvert Height : %s\n'% self.height
message += 'Culvert Width : %s\n'% self.width
message += 'Batter Slope 1 : %s\n'% self.z1
message += 'Batter Slope 2 : %s\n'% self.z2
message += 'Culvert Blockage: %s\n'% self.blockage
message += 'No. of barrels: %s\n'% self.barrels

#print "Structure Myids ",self.myid, self.label

Expand Down
79 changes: 69 additions & 10 deletions anuga/parallel/parallel_weir_orifice_trapezoid_operator.py
View file
Original file line number Diff line number Diff line change
@@ -1,5 +1,6 @@
import anuga
import math
import numpy

from anuga.structures.weir_orifice_trapezoid_operator import weir_orifice_trapezoid_function

Expand All @@ -25,14 +26,16 @@ def __init__(self,
z1=None,
z2=None,
blockage=0.0,
barrels=1.0,
end_points=None,
exchange_lines=None,
enquiry_points=None,
invert_elevations=None,
#culvert_slope=None,
apron=0.1,
manning=0.013,
enquiry_gap=0.0,
smoothing_timescale=0.0,
smoothing_timescale=0.0,
use_momentum_jet=True,
use_velocity_head=True,
description=None,
Expand All @@ -54,9 +57,11 @@ def __init__(self,
invert_elevations=invert_elevations,
width=width,
height=height,
blockage=blockage,
blockage=blockage,
barrels=barrels,
z1=z1,
z2=z2,
#culvert_slope=culvert_slope,
diameter= None,
apron=apron,
manning=manning,
Expand Down Expand Up @@ -92,6 +97,10 @@ def __init__(self,
self.culvert_length = self.get_culvert_length()
self.culvert_width = self.get_culvert_width()
self.culvert_height = self.get_culvert_height()
self.culvert_blockage = self.get_culvert_blockage()
self.culvert_barrels = self.get_culvert_barrels()

#self.culvert_slope = self.get_culvert_slope()

self.culvert_z1 = self.get_culvert_z1()
self.culvert_z2 = self.get_culvert_z2()
Expand All @@ -108,8 +117,8 @@ def __init__(self,

self.case = 'N/A'



self.domain=domain
# May/June 2014 -- allow 'smoothing ' of driving_energy, delta total energy, and outflow_enq_depth
self.smoothing_timescale=0.
self.smooth_delta_total_energy=0.
Expand All @@ -131,12 +140,18 @@ def __init__(self,


def parallel_safe(self):
"""
Set that operator is parallel safe
"""

return True



def discharge_routine(self):
"""
Get info from inlets and then call sequential function
"""

import pypar

Expand Down Expand Up @@ -181,19 +196,39 @@ def discharge_routine(self):
self.outflow_index = 1
# master proc orders reversal if applicable
if self.myid == self.master_proc:

# May/June 2014 -- change the driving forces gradually, with forward euler timestepping
#
forward_Euler_smooth=True
if(forward_Euler_smooth):
# To avoid 'overshoot' we ensure ts<1.
if(self.domain.timestep>0.):
ts=self.domain.timestep/max(self.domain.timestep, self.smoothing_timescale,1.0e-06)
else:
# This case is included in the serial version, which ensures the unit tests pass
# even when domain.timestep=0.0.
# Note though the discontinuous behaviour as domain.timestep-->0. from above
ts=1.0
self.smooth_delta_total_energy=self.smooth_delta_total_energy+\
ts*(self.delta_total_energy-self.smooth_delta_total_energy)
else:
# Use backward euler -- the 'sensible' ts limitation is different in this case
# ts --> Inf is reasonable and corresponds to the 'nosmoothing' case
ts=self.domain.timestep/max(self.smoothing_timescale, 1.0e-06)
self.smooth_delta_total_energy = (self.smooth_delta_total_energy+ts*(self.delta_total_energy))/(1.+ts)

# Reverse the inflow and outflow direction?
if self.delta_total_energy < 0:
if self.smooth_delta_total_energy < 0:
self.inflow_index = 1
self.outflow_index = 0

self.delta_total_energy = -self.delta_total_energy
#self.delta_total_energy = -self.delta_total_energy
self.delta_total_energy = -self.smooth_delta_total_energy

for i in self.procs:
if i == self.master_proc: continue
pypar.send(True, i)
else:
self.delta_total_energy = self.smooth_delta_total_energy
for i in self.procs:
if i == self.master_proc: continue
pypar.send(False, i)
Expand Down Expand Up @@ -261,22 +296,46 @@ def discharge_routine(self):
self.driving_energy = inflow_enq_specific_energy
else:
self.driving_energy = inflow_enq_depth


Q, barrel_velocity, outlet_culvert_depth, flow_area, case = \
weir_orifice_trapezoid_function(depth =self.culvert_height,
width =self.culvert_width,
z1 =self.culvert_z1,
z2 =self.culvert_z2,
flow_width =self.culvert_width,
length =self.culvert_length,
blockage =self.culvert_blockage,
barrels =self.culvert_barrels,
#culvert_slope =self.culvert_slope,
driving_energy =self.driving_energy,
delta_total_energy =self.delta_total_energy,
outlet_enquiry_depth=outflow_enq_depth,
sum_loss =self.sum_loss,
manning =self.manning)

################################################
# Smooth discharge. This can reduce oscillations
#
# NOTE: The sign of smooth_Q assumes that
# self.inflow_index=0 and self.outflow_index=1
# , whereas the sign of Q is always positive
Qsign=(self.outflow_index-self.inflow_index) # To adjust sign of Q
if(forward_Euler_smooth):
self.smooth_Q = self.smooth_Q +ts*(Q*Qsign-self.smooth_Q)
else:
# Try implicit euler method
self.smooth_Q = (self.smooth_Q+ts*(Q*Qsign))/(1.+ts)

if numpy.sign(self.smooth_Q)!=Qsign:
# The flow direction of the 'instantaneous Q' based on the
# 'smoothed delta_total_energy' is not the same as the
# direction of smooth_Q. To prevent 'jumping around', let's
# set Q to zero
Q=0.
else:
Q = min(abs(self.smooth_Q), Q) #abs(self.smooth_Q)
barrel_velocity=Q/flow_area
# END CODE BLOCK for DEPTH > Required depth for CULVERT Flow

else: # self.inflow.get_enquiry_depth() < 0.01:
Expand All @@ -291,7 +350,7 @@ def discharge_routine(self):
barrel_velocity = self.max_velocity
Q = flow_area * barrel_velocity



return Q, barrel_velocity, outlet_culvert_depth
else:
Expand Down
5 changes: 4 additions & 1 deletion anuga/structures/structure_operator.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -34,6 +34,7 @@ def __init__(self,
z2=None,
blockage=None,
barrels=None,
#culvert_slope=None,
apron=None,
manning=None,
enquiry_gap=None,
Expand Down Expand Up @@ -514,6 +515,8 @@ def statistics(self):
message += 'Culvert Width: %s\n'% self.width
message += 'Batter Slope 1: %s\n'% self.z1
message += 'Batter Slope 2: %s\n'% self.z2
message += 'Culvert Blockage: %s\n'% self.blockage
message += 'No. of barrels: %s\n'% self.barrels

message += '\n'

Expand Down Expand Up @@ -678,7 +681,7 @@ def get_culvert_blockage(self):
def get_culvert_barrels(self):

return self.barrels

def get_culvert_apron(self):

return self.apron
Expand Down

0 comments on commit 5859071

Please sign in to comment.