rhoPimpleFoamAvanzadoFlujo CompresibleOpenCFD·OpenFOAM v2606~55 minutos

Tubo de Choque de Sod con rhoPimpleFoam: ondas compresibles y discontinuidades

Caso unidimensional transitorio para observar una onda de choque, una discontinuidad de contacto y un abanico de rarefacción, con comparación cuantitativa frente a la solución de Riemann.

Qué aprenderás y validarás en este caso:

Configurar un caso compresible transitorio con p, U y T.
Inicializar dos estados termodinámicos mediante setFields.
Usar paso temporal adaptativo controlado por Courant.
Comparar presión, densidad y velocidad con una solución de referencia 1D.
Reconocer oscilaciones y difusión alrededor de discontinuidades.

Visualización de resultados CFD para Tubo de Choque de Sod con rhoPimpleFoam: ondas compresibles y discontinuidades

Estructura del caso OpenFOAM

Estructura del caso: rhoPimpleFoam

sodShockTube/
├── 0/
│   ├── p
│   ├── T
│   └── U
├── constant/
│   ├── thermophysicalProperties
│   └── turbulenceProperties
├── system/
│   ├── blockMeshDict
│   ├── controlDict
│   ├── fvSchemes
│   ├── fvSolution
│   └── setFieldsDict
├── Allrun
└── Allclean

Guía de ejecución paso a paso

Crear una malla 1D con resolución controlable

Construimos un dominio de longitud 1 m, una sola celda en y y z, y 400 celdas uniformes en x.

system/blockMeshDict(text)

1
2
3
4
5
6
7
8
9
10
scale 1;
vertices ((0 0 0) (1 0 0) (1 0.01 0) (0 0.01 0)
          (0 0 0.01) (1 0 0.01) (1 0.01 0.01) (0 0.01 0.01));
blocks (hex (0 1 2 3 4 5 6 7) (400 1 1) simpleGrading (1 1 1));
boundary
(
  left  { type patch; faces ((0 4 7 3)); }
  right { type patch; faces ((1 2 6 5)); }
  sides { type empty; faces ((0 1 5 4) (3 7 6 2) (0 3 2 1) (4 5 6 7)); }
);

Nota técnica: Una malla uniforme permite atribuir el ensanchamiento del choque al esquema y la resolución. Repetiremos con 200 y 800 celdas.

Definir los estados izquierdo y derecho

El estado base es p=0.1 Pa y T=0.8 K; setFields sobrescribe la mitad izquierda con p=1 Pa y T=1 K. Los valores están normalizados y se usa gas ideal.

system/setFieldsDict(text)

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
defaultFieldValues
(
    volScalarFieldValue p 0.1
    volScalarFieldValue T 0.8
    volVectorFieldValue U (0 0 0)
);
regions
(
    boxToCell
    {
        box (0 0 0) (0.5 0.01 0.01);
        fieldValues
        (
            volScalarFieldValue p 1.0
            volScalarFieldValue T 1.0
        );
    }
);

Nota técnica: Con R=1.25 en unidades normalizadas, ambos estados producen las densidades clásicas rho_L=1 y rho_R=0.1.

Configurar gas ideal y régimen laminar

Elegimos un gas caloricamente perfecto con gamma=1.4 y desactivamos turbulencia: el problema es un test inviscido de propagación de ondas.

constant/thermophysicalProperties(text)

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
thermoType
{
    type            hePsiThermo;
    mixture         pureMixture;
    transport       const;
    thermo          hConst;
    equationOfState perfectGas;
    specie          specie;
    energy          sensibleEnthalpy;
}
mixture
{
    specie { molWeight 6.6512; }
    thermodynamics { Cp 3.5; Hf 0; }
    transport { mu 1e-08; Pr 0.72; }
}

Nota técnica: La viscosidad pequeña reduce efectos físicos difusivos; el propósito es evaluar captura de discontinuidades, no modelar un gas real.

Controlar Courant y usar esquemas acotados

Activamos ajuste automático de deltaT con maxCo=0.4 y convección upwind robusta como referencia inicial.

system/controlDict(text)

1
2
3
4
5
6
7
8
application     rhoPimpleFoam;
endTime         0.2;
deltaT          1e-5;
adjustTimeStep  yes;
maxCo           0.4;
maxDeltaT       0.002;
writeControl    adjustableRunTime;
writeInterval   0.02;

Nota técnica: En flujo compresible la velocidad característica incluye U±a. Un deltaT basado solo en U subestima la restricción al inicio.

Muestrear el eje y comparar con la solución de Riemann

Extraemos p, U, T y rho sobre una línea. A t=0.2 deben distinguirse rarefacción a la izquierda, contacto central y choque a la derecha.

system/controlDict(text)

1
2
3
4
5
6
7
8
9
10
11
12
functions
{
    perfil
    {
        type sets;
        libs (sampling);
        writeControl writeTime;
        setFormat raw;
        fields (p U T rho);
        sets (eje { type uniform; axis x; start (0.001 0.005 0.005); end (0.999 0.005 0.005); nPoints 500; });
    }
}

Nota técnica: Compara posición de ondas y valores de meseta, no el valor exacto en una discontinuidad, que depende de resolución e interpolación.

Diccionarios y Archivos Completos del Caso

Explorador de Archivos Completos del Caso

Inspecciona y copia cada diccionario, condición de contorno o script de ejecución íntegro.

Archivos del Caso (12)

tubo-choque-sod-rhopimplefoam/system/blockMeshDict

10 líneas

system/blockMeshDict(text)

1
2
3
4
5
6
7
8
9
10
FoamFile { version 2.0; format ascii; class dictionary; object blockMeshDict; }
scale 1;
vertices ((0 0 0) (1 0 0) (1 0.01 0) (0 0.01 0) (0 0 0.01) (1 0 0.01) (1 0.01 0.01) (0 0.01 0.01));
blocks (hex (0 1 2 3 4 5 6 7) (400 1 1) simpleGrading (1 1 1));
boundary
(
 left { type patch; faces ((0 4 7 3)); }
 right { type patch; faces ((1 2 6 5)); }
 sides { type empty; faces ((0 1 5 4) (3 7 6 2) (0 3 2 1) (4 5 6 7)); }
);

Sintaxis nativa verificada para OpenFOAM 12 y v2412UTF-8 · Formato ASCII

Validación Cuantitativa vs Datos Experimentales / Canónicos

Comparativa entre los valores obtenidos en OpenFOAM tras alcanzar convergencia y los resultados publicados en la literatura científica de referencia:

Parámetro / Métrica	Dato Experimental / Referencia	Resultado OpenFOAM	Desviación Relativa
Presión detrás del choque	0.303 (solución exacta)	0.30–0.31	< 2%
Velocidad de la región estrella	0.927	0.91–0.94	< 2%
Posición del choque a t=0.2	x≈0.850	x≈0.85	< 1 celda fina

Guía de Análisis y Postprocesado en ParaView

Representa p, rho y U frente a x para 200, 400 y 800 celdas. Comprueba que la posición de las ondas converge y que el espesor numérico de la discontinuidad disminuye sin introducir sobreoscilaciones.

Problemas frecuentes y resolución de errores

⚠️ deltaT cae a valores extremadamente pequeños

Causa probable: Oscilaciones de presión o temperatura producen una velocidad del sonido no física.

Solución: Vuelve a upwind, revisa estados y thermo, limita maxCo a 0.2 durante el arranque y localiza el primer valor negativo.

⚠️ La densidad no aparece en sampling

Causa probable: rho puede no escribirse como campo según configuración.

Solución: Añade un function object de expresión/campo o reconstruye rho=p/(R*T) durante el postproceso.

Variaciones sugeridas del ejercicio

Comparar upwind con limitedLinear y medir espesor del choque.
Repetir con 200, 400 y 800 celdas.
Cambiar la relación de presiones para estudiar la intensidad del choque.

Referencias bibliográficas y de validación

OpenCFD documentation: rhoPimpleFoam solver overview.
Sod, G. A. (1978). A survey of several finite difference methods for systems of nonlinear hyperbolic conservation laws.

Volver al catálogo de tutoriales ¿Dudas sobre este caso? Pregunta en el foro