Cp (Pressure Coefficient)

Computation of pressure coefficient


Input Parameters

The following keys can be set:

  • base – (type = antares.Base ) – A base containing:

    • the mesh coordinates x, y, and z

    • the solution

    • ‘hb_computation’ as an Base.attrs (if HB/TSM type)

  • coordinates – (type = list(str) ) – the name of the variables that defines the mesh

  • vectors – (default = [], type = tuple/list of tuples of variables ) – if the base contains vectors, these must be rotated, so put them here. It is assumed that these are in cartesian coordinates

  • rho_inf – (default = in_attr, type = float, can use in_attr = yes ) – The infinite density corresponding to the reference state

  • v_inf – (default = in_attr, type = float, can use in_attr = yes ) – The infinite axial velocity corresponding to the reference state

  • p_inf – (default = in_attr, type = float, can use in_attr = yes ) – The infinite static pressure corresponding to the reference state

  • family_name – (type = str ) – The name of the family from which the percent will be computed and on which Cp is computed

  • percent – (default = None, type = float or None ) – The percentage relative to the family

  • position – (default = None, type = float or None ) – The absolute position value relative to the family where the cut must be made

  • form – (default = 1, type = int in [1,2,3] ) – the definition of Cp (see below)

Main functions

class antares.treatment.turbomachine.TreatmentCp.TreatmentCp
execute()

Execute the treatment.

Compute the pressure coefficient at a given radius percent of the given blade. Thee formulae for Cp are proposed:

form 1: \(\displaystyle Cp_1 = - \frac{p - p_{inf}}{\rho_{inf} n^2 D^2}\)

form 2: \(\displaystyle Cp_2 = 2 \frac{p - p_{inf}}{\rho_{inf} (v_{inf} n^2 r^2)}\)

form 3: \(\displaystyle Cp_3 = 2.0 \frac{p - p_{inf}} {\rho_{inf} v_{inf} ^ 2}\).

The mean and the harmonics of Cp can also be computed if duplication is enabled. Note that the amplitude of the harmonics are divided by the mean value.

Returns

Return type

antares.Base

Example

import os

if not os.path.isdir('OUTPUT'):
    os.makedirs('OUTPUT')

from antares import Reader, Treatment, Writer

#

# Data can be downloaded from
# https://cerfacs.fr/antares/tutorial/application/application1/application1_tutorial_data.tgz

r = Reader('bin_tp')
r['filename'] = os.path.join('..', 'data', 'ROTOR37', 'ELSA_CASE', 'MESH',
                             'mesh_<zone>.dat')
r['zone_prefix'] = 'Block'
r['topology_file'] = os.path.join('..', 'data', 'ROTOR37', 'ELSA_CASE',
                                  'script_topo.py')
r['shared'] = True
base = r.read()
print(base.families)

r = Reader('bin_tp')
r['base'] = base
r['filename'] = os.path.join('..', 'data', 'ROTOR37', 'ELSA_CASE', 'FLOW',
                             'flow_<zone>.dat')
r['zone_prefix'] = 'Block'
r['location'] = 'cell'
r.read()

base.set_computer_model('internal')

# Needed for turbomachinery dedicated treatments
base.cell_to_node()
base = base.get_location('node')
print(base.families)

base.compute('psta')
base.compute('Pi')
base.compute('theta')
P0_INF = 1.9
base.compute('MachIs = (((%f/psta)**((gamma-1)/gamma)-1.) * (2./(gamma-1.))  )**0.5' % P0_INF)

# Definition of the treatment
t = Treatment('Cp')
t['base'] = base
t['family_name'] = 'BLADE'
t['coordinates'] = ['x', 'y', 'z']
t['rho_inf'] = 0.873
t['p_inf'] = 0.59
t['v_inf'] = 1.5
t['form'] = 3

# Cp
res_dir = os.path.join('OUTPUT', 'CP')
if not os.path.isdir(res_dir):
    os.makedirs(res_dir)

writer = Writer('column')

for loc in [0.25, 0.5, 0.75, 0.9]:  # radius in percent
    t['percent'] = loc
    Cp_blade = t.execute()

    writer['filename'] = os.path.join(res_dir, 'Cp_%s.dat' % (loc))
    writer['base'] = Cp_blade
    writer.dump()