Skip to content

Commit

Permalink
Added unit tests for conical grains
Browse files Browse the repository at this point in the history
  • Loading branch information
@reilleya
reilleya committed Oct 2, 2021
1 parent 7c1b1d1 commit 9e4d6fd
Show file tree
Hide file tree
Showing 3 changed files with 206 additions and 49 deletions.
14 changes: 7 additions & 7 deletions motorlib/geometry.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -27,21 +27,21 @@ def cylinderVolume(dia, height):
"""Returns the volume of a cylinder with diameter 'dia' and height 'height'"""
return height * circleArea(dia)

def frustrumLateralSurfaceArea(diameterA, diameterB, length):
"""Returns the surface area of a frustrum (truncated cone) with end diameters A and B and length 'length'"""
def frustumLateralSurfaceArea(diameterA, diameterB, length):
"""Returns the surface area of a frustum (truncated cone) with end diameters A and B and length 'length'"""
radiusA = diameterA / 2
radiusB = diameterB / 2
return math.pi * (radiusA + radiusB) * (abs(radiusA - radiusB) ** 2 + length ** 2) ** 0.5

def frustrumVolume(diameterA, diameterB, length):
"""Returns the volume of a frustrum (truncated cone) with end diameters A and B and length 'length'"""
def frustumVolume(diameterA, diameterB, length):
"""Returns the volume of a frustum (truncated cone) with end diameters A and B and length 'length'"""
radiusA = diameterA / 2
radiusB = diameterB / 2
return math.pi * (length / 3) * (radiusA ** 2 + radiusA * radiusB + radiusB ** 2)

def splitFrustrum(diameterA, diameterB, length, splitPosition):
"""Takes in info about a frustrum (truncated cone) and a position measured from the "diameterA" and returns
a tuple of frustrums representing the two halves of the original frustrum if it were split on the plane at
def splitFrustum(diameterA, diameterB, length, splitPosition):
"""Takes in info about a frustum (truncated cone) and a position measured from the "diameterA" and returns
a tuple of frustums representing the two halves of the original frustum if it were split on the plane at
distance "position" from the face with diameter "diameterA"
"""
splitDiameter = diameterA + (diameterB - diameterA) * (splitPosition / length)
Expand Down
54 changes: 28 additions & 26 deletions motorlib/grains/conical.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -23,8 +23,8 @@ def isCoreInverted(self):
"""A simple helper that returns 'true' if the core's foward diameter is larger than its aft diameter"""
return self.props['forwardCoreDiameter'].getValue() > self.props['aftCoreDiameter'].getValue()

def getFrustrumInfo(self, regDist):
"""Returns the dimensions of the grain's core at a given regression depth. The core is always a frustrum and is
def getFrustumInfo(self, regDist):
"""Returns the dimensions of the grain's core at a given regression depth. The core is always a frustum and is
returned as the aft diameter, forward diameter, and length"""
grainDiameter = self.props['diameter'].getValue()
aftDiameter = self.props['aftCoreDiameter'].getValue()
Expand Down Expand Up @@ -57,11 +57,11 @@ def getFrustrumInfo(self, regDist):
# the diameter of the large end of the core is clamped at the grain diameter and the length is changed to keep
# the angle constant, which accounts for the regression of the grain at the major end.
if regCoreMajorDiameter >= grainDiameter:
majorFrustrumDiameter = grainDiameter
minorFrustrumDiameter = regCoreMinorDiameter
majorFrustumDiameter = grainDiameter
minorFrustumDiameter = regCoreMinorDiameter
# How far past the grain diameter the major end has grown
effectiveReg = (regCoreMajorDiameter - grainDiameter) / 2
# Reduce the length of the frustrum by the axial component of the regression vector
# Reduce the length of the frustum by the axial component of the regression vector
grainLength -= (effectiveReg / sin(angle))
# Account for the change in length due to face regression before the major end hit the casting tube
grainLength -= exposedFaces * (grainDiameter - coreMajorDiameter) / 2
Expand All @@ -71,19 +71,19 @@ def getFrustrumInfo(self, regDist):
# If the large end of the core hasn't reached the casting tube, we know that the small end hasn't either. In
# this case we just return the current core diameters, and a length calculated from the inhibitor configuration
else:
majorFrustrumDiameter = regCoreMajorDiameter
minorFrustrumDiameter = regCoreMinorDiameter
majorFrustumDiameter = regCoreMajorDiameter
minorFrustumDiameter = regCoreMinorDiameter
grainLength -= exposedFaces * regDist

if self.isCoreInverted():
return minorFrustrumDiameter, majorFrustrumDiameter, grainLength
return minorFrustumDiameter, majorFrustumDiameter, grainLength

return majorFrustrumDiameter, minorFrustrumDiameter, grainLength
return majorFrustumDiameter, minorFrustumDiameter, grainLength

def getSurfaceAreaAtRegression(self, regDist):
"""Returns the surface area of the grain after it has regressed a linear distance of 'regDist'"""
aftDiameter, forwardDiameter, length = self.getFrustrumInfo(regDist)
surfaceArea = geometry.frustrumLateralSurfaceArea(aftDiameter, forwardDiameter, length)
aftDiameter, forwardDiameter, length = self.getFrustumInfo(regDist)
surfaceArea = geometry.frustumLateralSurfaceArea(aftDiameter, forwardDiameter, length)

fullFaceArea = geometry.circleArea(self.props['diameter'].getValue())
if self.props['inhibitedEnds'].getValue() in ['Neither', 'Bottom']:
Expand All @@ -95,15 +95,15 @@ def getSurfaceAreaAtRegression(self, regDist):

def getVolumeAtRegression(self, regDist):
"""Returns the volume of propellant in the grain after it has regressed a linear distance 'regDist'"""
aftDiameter, forwardDiameter, length = self.getFrustrumInfo(regDist)
frustrumVolume = geometry.frustrumVolume(aftDiameter, forwardDiameter, length)
aftDiameter, forwardDiameter, length = self.getFrustumInfo(regDist)
frustumVolume = geometry.frustumVolume(aftDiameter, forwardDiameter, length)
outerVolume = geometry.cylinderVolume(self.props['diameter'].getValue(), length)

return outerVolume - frustrumVolume
return outerVolume - frustumVolume

def getWebLeft(self, regDist):
"""Returns the shortest distance the grain has to regress to burn out"""
aftDiameter, forwardDiameter, length = self.getFrustrumInfo(regDist)
aftDiameter, forwardDiameter, length = self.getFrustumInfo(regDist)

return (self.props['diameter'].getValue() - min(aftDiameter, forwardDiameter)) / 2

Expand All @@ -113,23 +113,25 @@ def getMassFlow(self, massIn, dTime, regDist, dRegDist, position, density):
position along the grain measured from the head end, and the density of the propellant."""

# For now these grains are only allowed with inhibited faces, so we can ignore a lot of messy logic
unsteppedFrustrum = self.getFrustrumInfo(regDist)
steppedFrustrum = self.getFrustrumInfo(regDist + dRegDist)
unsteppedFrustum = self.getFrustumInfo(regDist)
steppedFrustum = self.getFrustumInfo(regDist + dRegDist)

unsteppedFrustrum = (unsteppedFrustrum[1], unsteppedFrustrum[0], unsteppedFrustrum[2])
steppedFrustrum = (steppedFrustrum[1], steppedFrustrum[0], steppedFrustrum[2])
"""These have to be reordered, because getFrustumInfo returns (aft, forward, length), but splitFrustum wants
the position from the forward face"""
unsteppedFrustum = (unsteppedFrustum[1], unsteppedFrustum[0], unsteppedFrustum[2])
steppedFrustum = (steppedFrustum[1], steppedFrustum[0], steppedFrustum[2])

# Note that this assumes the forward end of the grain is still at postition 0 - inhibited
unsteppedPartialFrustrum, _ = geometry.splitFrustrum(*unsteppedFrustrum, position)
steppedPartialFrustrum, _ = geometry.splitFrustrum(*steppedFrustrum, position)
unsteppedPartialFrustum, _ = geometry.splitFrustum(*unsteppedFrustum, position)
steppedPartialFrustum, _ = geometry.splitFrustum(*steppedFrustum, position)

unsteppedVolume = geometry.frustrumVolume(*unsteppedPartialFrustrum)
steppedVolume = geometry.frustrumVolume(*steppedPartialFrustrum)
unsteppedVolume = geometry.frustumVolume(*unsteppedPartialFrustum)
steppedVolume = geometry.frustumVolume(*steppedPartialFrustum)

massFlow = (steppedVolume - unsteppedVolume) * density / dTime
massFlow += massIn

return massFlow, steppedPartialFrustrum[1]
return massFlow, steppedPartialFrustum[1]

def getMassFlux(self, massIn, dTime, regDist, dRegDist, position, density):
"""Returns the mass flux at a point along the grain. Takes in the mass flow into the grain, a timestep, the
Expand All @@ -153,7 +155,7 @@ def getEndPositions(self, regDist):
"""Returns the positions of the grain ends relative to the original (unburned) grain top. Returns a tuple like
(forward, aft)"""
originalLength = self.props['length'].getValue()
aftCoreDiameter, forwardCoreDiameter, currentLength = self.getFrustrumInfo(regDist)
aftCoreDiameter, forwardCoreDiameter, currentLength = self.getFrustumInfo(regDist)
inhibitedEnds = self.props['inhibitedEnds'].getValue()

if self.isCoreInverted():
Expand All @@ -177,7 +179,7 @@ def getEndPositions(self, regDist):

def getPortArea(self, regDist):
"""Returns the area of the grain's port when it has regressed a distance of 'regDist'"""
aftCoreDiameter, _, _ = self.getFrustrumInfo(regDist)
aftCoreDiameter, _, _ = self.getFrustumInfo(regDist)

return geometry.circleArea(aftCoreDiameter)

Expand Down
187 changes: 171 additions & 16 deletions test/test.py
View file
Original file line number Diff line number Diff line change
@@ -1,5 +1,10 @@
import unittest


#######################################################################################################################
### Motor Tests ###
#######################################################################################################################

import motorlib.motor
import motorlib.grains
import motorlib.propellant
Expand All @@ -10,10 +15,10 @@ def test_calcKN(self):
tc = motorlib.motor.MotorConfig()

bg = motorlib.grains.BatesGrain()
bg.setProperties({'diameter':0.083058,
'length':0.1397,
'coreDiameter':0.05,
'inhibitedEnds':'Neither'
bg.setProperties({'diameter': 0.083058,
'length': 0.1397,
'coreDiameter': 0.05,
'inhibitedEnds': 'Neither'
})

tm.grains.append(bg)
Expand All @@ -24,17 +29,17 @@ def test_calcKN(self):
self.assertAlmostEqual(tm.calcKN([0.0025], 0), 183, 0)
self.assertAlmostEqual(tm.calcKN([0.005], 0), 185, 0)


def test_calcPressure(self):
tm = motorlib.motor.Motor()
tc = motorlib.motor.MotorConfig()

bg = motorlib.grains.BatesGrain()
bg.setProperties({'diameter':0.083058,
'length':0.1397,
'coreDiameter':0.05,
'inhibitedEnds':'Neither'
})
bg.setProperties({
'diameter': 0.083058,
'length': 0.1397,
'coreDiameter': 0.05,
'inhibitedEnds': 'Neither'
})

tm.grains.append(bg)
bg.simulationSetup(tc)
Expand All @@ -43,19 +48,24 @@ def test_calcPressure(self):
tm.propellant = motorlib.propellant.Propellant()
tm.propellant.setProperties({
'name': 'KNSU',
'density': 1890,
'tabs':[
'density': 1890,
'tabs': [
{
'a': 0.000101,
'n': 0.319,
't': 1720,
'm': 41.98,
'a': 0.000101,
'n': 0.319,
't': 1720,
'm': 41.98,
'k': 1.133
}
]
})
self.assertAlmostEqual(tm.calcIdealPressure([0], 0), 4050196, 0)


#######################################################################################################################
### Geometry Tests ###
#######################################################################################################################

import motorlib.geometry
class TestGeometryMethods(unittest.TestCase):
def test_circleArea(self):
Expand All @@ -76,6 +86,24 @@ def test_cylinderArea(self):
def test_cylinderVolume(self):
self.assertAlmostEqual(motorlib.geometry.cylinderVolume(0.5, 2), 0.39269908)

def test_frustumLateralSurfaceArea(self):
self.assertAlmostEqual(motorlib.geometry.frustumLateralSurfaceArea(2, 3, 5), 39.46576927)

def test_frustumVolume(self):
# Cone case
self.assertAlmostEqual(motorlib.geometry.frustumVolume(0, 10, 10), 261.79938779)
# Frustum case
self.assertAlmostEqual(motorlib.geometry.frustumVolume(10, 30, 50), 17016.96020694)

def test_splitFrustum(self):
# Simple case
self.assertAlmostEqual(motorlib.geometry.splitFrustum(1, 2, 4, 2), ((1, 1.5, 2), (1.5, 2, 2)))
# Inverted case
self.assertAlmostEqual(motorlib.geometry.splitFrustum(2, 1, 4, 2), ((2, 1.5, 2), (1.5, 1, 2)))
# Make sure that the connected ends of the frustums line up
upper, lower = motorlib.geometry.splitFrustum(1, 3, 3, 1)
self.assertEqual(upper[1], lower[0])

def test_dist(self):
self.assertEqual(motorlib.geometry.dist((5, 5), (5, 5)), 0)
self.assertEqual(motorlib.geometry.dist((5, 5), (6, 5)), 1)
Expand Down Expand Up @@ -115,6 +143,11 @@ def test_getExitPressure(self):
self.assertAlmostEqual(nozzle.getExitPressure(1.2, 5e6), 72087.22454540983)
self.assertAlmostEqual(nozzle.getExitPressure(1.2, 6e6), 86504.66945449157)


#######################################################################################################################
### Propellant Tests ###
#######################################################################################################################

import motorlib.propellant
class TestPropellantMethods(unittest.TestCase):
def test_proper_propellant_ranges(self):
Expand Down Expand Up @@ -245,4 +278,126 @@ def test_get_combustion_properties_out_of_range(self):
self.assertEqual(testProp.getCombustionProperties(6.9e5), (1.467e-05, 0.382, 1.25, 3500, 23.67))
self.assertEqual(testProp.getCombustionProperties(8e10), (1e-05, 0.3, 1.25, 3500, 23.67))


#######################################################################################################################
### Conical Grain Tests ###
#######################################################################################################################

import motorlib.grains
class ConicalGrainMethods(unittest.TestCase):

def test_isCoreInverted(self):
inverted = motorlib.grains.ConicalGrain()
inverted.setProperties({
'length': 0.1,
'diameter': 0.01,
'forwardCoreDiameter': 0.0025,
'aftCoreDiameter': 0.002,
})
regular = motorlib.grains.ConicalGrain()
regular.setProperties({
'length': 0.1,
'diameter': 0.01,
'forwardCoreDiameter': 0.003,
'aftCoreDiameter': 0.004,
})

self.assertEqual(inverted.isCoreInverted(), True)
self.assertEqual(regular.isCoreInverted(), False)

def test_getFrustumInfo(self):
properties = {
'length': 0.1,
'diameter': 0.01,
'forwardCoreDiameter': 0.0025,
'aftCoreDiameter': 0.002,
'inhibitedEnds': 'Both'
}

testGrain = motorlib.grains.ConicalGrain()
testGrain.setProperties(properties)

unregressed = testGrain.getFrustumInfo(0)
self.assertAlmostEqual(unregressed[0], properties['aftCoreDiameter'])
self.assertAlmostEqual(unregressed[1], properties['forwardCoreDiameter'])
self.assertAlmostEqual(unregressed[2], properties['length'])

beforeHittingWall = testGrain.getFrustumInfo(0.001)
self.assertAlmostEqual(beforeHittingWall[0], 0.003999993750029297)
self.assertAlmostEqual(beforeHittingWall[1], 0.004499993750029296)
self.assertAlmostEqual(beforeHittingWall[2], properties['length']) # Length hasn't changed yet

hitWall = testGrain.getFrustumInfo(0.0038)
self.assertAlmostEqual(hitWall[0], 0.009599976250111327)
self.assertAlmostEqual(hitWall[1], properties['diameter']) # This end has burned all the way to the wall
self.assertAlmostEqual(hitWall[2], 0.08000468749267584)

def test_getSurfaceAreaAtRegression(self):
properties = {
'length': 0.1,
'diameter': 0.01,
'forwardCoreDiameter': 0.0025,
'aftCoreDiameter': 0.002,
'inhibitedEnds': 'Both'
}

forwardFaceArea = 7.36310778e-05
aftFaceArea = 7.53982236e-05
lateralArea = 0.00070686055598659

testGrain = motorlib.grains.ConicalGrain()
testGrain.setProperties(properties)

self.assertAlmostEqual(testGrain.getSurfaceAreaAtRegression(0), lateralArea)
self.assertAlmostEqual(testGrain.getSurfaceAreaAtRegression(0.001), 0.0013351790867045452)

# For when uninibited conical grains work:
"""testGrain.setProperty('inhibitedEnds', 'Top')
self.assertAlmostEqual(testGrain.getSurfaceAreaAtRegression(0), lateralArea + aftFaceArea)
testGrain.setProperty('inhibitedEnds', 'Bottom')
self.assertAlmostEqual(testGrain.getSurfaceAreaAtRegression(0), lateralArea + forwardFaceArea)
testGrain.setProperty('inhibitedEnds', 'Neither')
self.assertAlmostEqual(testGrain.getSurfaceAreaAtRegression(0), lateralArea + forwardFaceArea + aftFaceArea)"""

def test_getVolumeAtRegression(self):
properties = {
'length': 0.1,
'diameter': 0.01,
'forwardCoreDiameter': 0.0025,
'aftCoreDiameter': 0.002,
'inhibitedEnds': 'Both'
}

testGrain = motorlib.grains.ConicalGrain()
testGrain.setProperties(properties)

self.assertAlmostEqual(testGrain.getVolumeAtRegression(0), 7.454737567580781e-06)
self.assertAlmostEqual(testGrain.getVolumeAtRegression(0.001), 6.433724127569215e-06)
self.assertAlmostEqual(testGrain.getVolumeAtRegression(0.0038), 2.480054353678591e-07)

def test_getWebLeft(self):
properties = {
'length': 0.1,
'diameter': 0.01,
'forwardCoreDiameter': 0.0025,
'aftCoreDiameter': 0.002,
'inhibitedEnds': 'Both'
}

testGrain = motorlib.grains.ConicalGrain()
testGrain.setProperties(properties)

self.assertAlmostEqual(testGrain.getWebLeft(0), 0.004)
self.assertAlmostEqual(testGrain.getWebLeft(0.001), 0.003)
self.assertAlmostEqual(testGrain.getWebLeft(0.0038), 0.0002)

testGrain.setProperty('forwardCoreDiameter', 0.002)
testGrain.setProperty('aftCoreDiameter', 0.0025)
self.assertAlmostEqual(testGrain.getWebLeft(0), 0.004)
self.assertAlmostEqual(testGrain.getWebLeft(0.001), 0.003)
self.assertAlmostEqual(testGrain.getWebLeft(0.0038), 0.0002)


unittest.main()

0 comments on commit 9e4d6fd

Please sign in to comment.