In [ ]:

    

%matplotlib inline
import numpy as np
import matplotlib
import matplotlib.pyplot as plt
vs  = np.linspace(0.0,1500,100.0);
TDs = np.linspace(0.25,4.0,0.25)
D0 = 100.0e-3

angles = np.linspace(0.0,np.pi*0.5*0.9, 100)
cosAs  = np.cos(angles) 
lds    = np.linspace(1.0,1.0, 100)

rad2deg = 180.0/np.pi


    

In [ ]:

    

def area( caliber ):
    return np.pi*(caliber*0.5)**2

def projectileMass( L, caliber=100e-3, density=7.8e+3 ):
    return area(caliber)*L*density

def kineticEnergy( v, mass ):
    return 0.5*mass*v**2

def flatIncidence( v, thickNor, caliber=100e-3, mass=None, density=7.8e+3, strength=2e+9, cShear=0.75, cDisp=1.0, debug=False ):
    if mass is None:
        mass = projectileMass( caliber*2, caliber=caliber, density=7.8e+3 )
        print "mass=", mass
    thick       =  thickNor*caliber
    #print "thick=", thick
    Ek          = kineticEnergy( v, mass )
    S           = area(caliber)
    Sshear      = np.pi*caliber*thick*cShear
    #print S,Sshear
    S           = np.minimum(S,Sshear) 
    Estatic        = thick*S*strength
    #print Ek,Estatic
    displaced_mass = thick*S*density
    Edynamic       = cDisp * kineticEnergy( v, displaced_mass )
    Eout       = np.maximum( Ek - Estatic - Edynamic , 0 )
    if debug:
        plt.plot(vs,Ek*1e-6,              label='Ek')
        plt.plot(vs,Edynamic*1e-6,        label='Edyn')
        plt.plot(vs,(Estatic+Ek*0)*1e-6,  label='Estat')
        plt.plot(vs,Eout*1e-6,'k-',lw=2,  label='Eout')
    return Eout

projectile_mass = projectileMass( 2*D0, caliber=D0, density=7.8e+3 )
print "projectile mass %g [kg]" %projectile_mass 

Eout = flatIncidence( vs, 1.0, caliber=D0, mass = projectile_mass, debug=True )
#plt.plot(vs,Eout*1e-6)

#flatIncidence( vs, projectile_mass, D0, caliber=D0, density=7.8e+3, strength=2e+9, cShear=0.75, cDisp=1.0 )

plt.xlabel("v [m/s]"); plt.ylabel( "E [MJ]" ); plt.legend(loc=2); plt.grid();


    

In [ ]:

    

def obiqueIncidence_1( v, thickNor, cosA ):
    Eout = flatIncidence( v, thickNor/cosA  )
    return Eout

#print cosAs
Eout = obiqueIncidence_1( 800.0, 1.0, cosAs )
plt.plot(angles*rad2deg, Eout*1e-6)

plt.xlabel("angle [deg]"); plt.ylabel( "Eout [MJ]" ); plt.legend(loc=4); plt.grid();


    

In [ ]:

    

def ricochetAngle_1( v, ld, Yp=2e+9, rho_p=7.8e+3, rho_t=7.8e+3 ):
    '''
    A simple estimate of the minimum target obliquity required for the ricochet of a high speed long rod projectile
    1979 J. Phys. D: Appl. Phys. 12 1825; http://iopscience.iop.org/0022-3727/12/11/011
    
    modified hydrodynamic equation:
    rho_p*(Vin-Vout)**2 + Yp = rho_t*Vout**2 + Yt
    Vcrit = sqrt( 2*(Rt-Yp)/rho_p )
    '''
    cDens  =  1.0 + np.sqrt(rho_p/rho_t)
    #cShape = (l**2+d**2)/(l*d)
    cShape = (ld**2 + 1)/(ld*1.0)
    print cDens, cShape
    rhs    = (2.0/3.0)*((rho_p * v**2)/Yp)*cShape*cDens
    tanA   = rhs**(1.0/3.0) 
    return np.arctan( tanA )

for ld in [1.0,2.0,3.0,4.0,5.0]:
    aCrit = ricochetAngle_1( vs, ld )
    plt.plot(vs,aCrit*rad2deg, label=("l/d=%1.1f" %ld) )

plt.xlabel("v [m/s]"); plt.ylabel( "aCrit [deg]" ); plt.legend(loc=4); plt.grid();


    

In [ ]:

    

def surfaceDeflection( v, angle, caliber=100e-3, mass=None, strength=2.0e+9 ):
    '''
    before projectile dent into the target it is accelerated along the surface
    '''
    tgA = np.tan(angle)
    #print tgA
    if mass is None:
        mass = projectileMass( caliber*2, caliber=caliber, density=7.8e+3 )
        print "mass=", mass
    tdent = (caliber*tgA)/v
    S   = area(caliber)
    aT  = ( strength * S * tgA ) / mass
    vT  = tdent * aT
    #print vT
    #vII = np.sqrt( vT )
    angle_ = angle  + np.arcsin( vT/v ) 
    return angle_

for v in [500,750,1000,1250]:
    angles_ = surfaceDeflection( v, angles )
    plt.plot( angles*rad2deg, angles_*rad2deg, label=("%1.1f [km/s]" %(v*1e-3)) )

plt.plot(angles*rad2deg,angles*rad2deg,ls='--',c="#808080")

plt.xlabel("alpha [deg]"); plt.ylabel( "beta [deg]" ); plt.legend(loc=4); plt.grid(); #plt.axis('equal')


    

In [ ]: