1## libCEED: Navier-Stokes Example 2 3This page provides a description of the Navier-Stokes example for the libCEED library, based on PETSc. 4PETSc v3.17 or a development version of PETSc at commit 0e95d842 or later is required. 5 6The Navier-Stokes problem solves the compressible Navier-Stokes equations in three dimensions using an explicit time integration. 7The state variables are mass density, momentum density, and energy density. 8 9The main Navier-Stokes solver for libCEED is defined in [`navierstokes.c`](navierstokes.c) with different problem definitions according to the application of interest. 10 11Build by using: 12 13`make` 14 15and run with: 16 17``` 18./navierstokes -ceed [ceed] -problem [problem type] -degree [degree] 19``` 20 21## Runtime options 22 23% inclusion-fluids-marker 24 25The Navier-Stokes mini-app is controlled via command-line options. 26The following options are common among all problem types: 27 28:::{list-table} Common Runtime Options 29:header-rows: 1 30 31* - Option 32 - Description 33 - Default value 34 35* - `-ceed` 36 - CEED resource specifier 37 - `/cpu/self/opt/blocked` 38 39* - `-test` 40 - Run in test mode 41 - `false` 42 43* - `-compare_final_state_atol` 44 - Test absolute tolerance 45 - `1E-11` 46 47* - `-compare_final_state_filename` 48 - Test filename 49 - 50 51* - `-problem` 52 - Problem to solve (`advection`, `advection2d`, `density_current`, or `euler_vortex`) 53 - `density_current` 54 55* - `-implicit` 56 - Use implicit time integartor formulation 57 - 58 59* - `-degree` 60 - Polynomial degree of tensor product basis (must be >= 1) 61 - `1` 62 63* - `-q_extra` 64 - Number of extra quadrature points 65 - `0` 66 67* - `-viz_refine` 68 - Use regular refinement for visualization 69 - `0` 70 71* - `-output_freq` 72 - Frequency of output, in number of steps. `0` has no output, `-1` outputs final state only 73 - `10` 74 75* - `-output_dir` 76 - Output directory 77 - `.` 78 79* - `-output_add_stepnum2bin` 80 - Whether to add step numbers to output binary files 81 - `false` 82 83* - `-continue` 84 - Continue from previous solution (input is step number of previous solution) 85 - `0` 86 87* - `-continue_filename` 88 - Path to solution binary file from which to continue from 89 - `[output_dir]/ns-solution.bin` 90 91* - `-continue_time_filename` 92 - Path to time stamp binary file from which to continue from 93 - `[output_dir]/ns-time.bin` 94 95* - `-bc_wall` 96 - Use wall boundary conditions on this list of faces 97 - 98 99* - `-wall_comps` 100 - An array of constrained component numbers for wall BCs 101 - 102 103* - `-bc_slip_x` 104 - Use slip boundary conditions, for the x component, on this list of faces 105 - 106 107* - `-bc_slip_y` 108 - Use slip boundary conditions, for the y component, on this list of faces 109 - 110 111* - `-bc_slip_z` 112 - Use slip boundary conditions, for the z component, on this list of faces 113 - 114 115* - `-bc_inflow` 116 - Use inflow boundary conditions on this list of faces 117 - 118 119* - `-bc_outflow` 120 - Use outflow boundary conditions on this list of faces 121 - 122 123* - `-bc_freestream` 124 - Use freestream boundary conditions on this list of faces 125 - 126 127* - `-snes_view` 128 - View PETSc `SNES` nonlinear solver configuration 129 - 130 131* - `-log_view` 132 - View PETSc performance log 133 - 134 135* - `-help` 136 - View comprehensive information about run-time options 137 - 138::: 139 140For the case of a square/cubic mesh, the list of face indices to be used with `-bc_wall`, `bc_inflow`, `bc_outflow`, `bc_freestream` and/or `-bc_slip_x`, `-bc_slip_y`, and `-bc_slip_z` are: 141 142:::{list-table} 2D Face ID Labels 143:header-rows: 1 144* - PETSc Face Name 145 - Cartesian direction 146 - Face ID 147 148* - faceMarkerBottom 149 - -z 150 - 1 151 152* - faceMarkerRight 153 - +x 154 - 2 155 156* - faceMarkerTop 157 - +z 158 - 3 159 160* - faceMarkerLeft 161 - -x 162 - 4 163::: 164 165:::{list-table} 2D Face ID Labels 166:header-rows: 1 167* - PETSc Face Name 168 - Cartesian direction 169 - Face ID 170 171* - faceMarkerBottom 172 - -z 173 - 1 174 175* - faceMarkerTop 176 - +z 177 - 2 178 179* - faceMarkerFront 180 - -y 181 - 3 182 183* - faceMarkerBack 184 - +y 185 - 4 186 187* - faceMarkerRight 188 - +x 189 - 5 190 191* - faceMarkerLeft 192 - -x 193 - 6 194::: 195 196### Advection 197 198For testing purposes, there is a reduced mode for pure advection, which holds density $\rho$ and momentum density $\rho \bm u$ constant while advecting "total energy density" $E$. 199These are available in 2D and 3D. 200 201#### 2D advection 202 203For the 2D advection problem, the following additional command-line options are available: 204 205:::{list-table} Advection2D Runtime Options 206:header-rows: 1 207 208* - Option 209 - Description 210 - Default value 211 - Unit 212 213* - `-rc` 214 - Characteristic radius of thermal bubble 215 - `1000` 216 - `m` 217 218* - `-units_meter` 219 - 1 meter in scaled length units 220 - `1E-2` 221 - 222 223* - `-units_second` 224 - 1 second in scaled time units 225 - `1E-2` 226 - 227 228* - `-units_kilogram` 229 - 1 kilogram in scaled mass units 230 - `1E-6` 231 - 232 233* - `-strong_form` 234 - Strong (1) or weak/integrated by parts (0) residual 235 - `0` 236 - 237 238* - `-stab` 239 - Stabilization method (`none`, `su`, or `supg`) 240 - `none` 241 - 242 243* - `-CtauS` 244 - Scale coefficient for stabilization tau (nondimensional) 245 - `0` 246 - 247 248* - `-wind_type` 249 - Wind type in Advection (`rotation` or `translation`) 250 - `rotation` 251 - 252 253* - `-wind_translation` 254 - Constant wind vector when `-wind_type translation` 255 - `1,0,0` 256 - 257 258* - `-E_wind` 259 - Total energy of inflow wind when `-wind_type translation` 260 - `1E6` 261 - `J` 262::: 263 264An example of the `rotation` mode can be run with: 265 266``` 267./navierstokes -problem advection2d -dm_plex_box_faces 20,20 -dm_plex_box_lower 0,0 -dm_plex_box_upper 1000,1000 -bc_wall 1,2,3,4 -wall_comps 4 -wind_type rotation -implicit -stab supg 268``` 269 270and the `translation` mode with: 271 272``` 273./navierstokes -problem advection2d -dm_plex_box_faces 20,20 -dm_plex_box_lower 0,0 -dm_plex_box_upper 1000,1000 -units_meter 1e-4 -wind_type translation -wind_translation 1,-.5 -bc_inflow 1,2,3,4 274``` 275Note the lengths in `-dm_plex_box_upper` are given in meters, and will be nondimensionalized according to `-units_meter`. 276 277#### 3D advection 278 279For the 3D advection problem, the following additional command-line options are available: 280 281:::{list-table} Advection3D Runtime Options 282:header-rows: 1 283 284* - Option 285 - Description 286 - Default value 287 - Unit 288 289* - `-rc` 290 - Characteristic radius of thermal bubble 291 - `1000` 292 - `m` 293 294* - `-units_meter` 295 - 1 meter in scaled length units 296 - `1E-2` 297 - 298 299* - `-units_second` 300 - 1 second in scaled time units 301 - `1E-2` 302 - 303 304* - `-units_kilogram` 305 - 1 kilogram in scaled mass units 306 - `1E-6` 307 - 308 309* - `-strong_form` 310 - Strong (1) or weak/integrated by parts (0) residual 311 - `0` 312 - 313 314* - `-stab` 315 - Stabilization method (`none`, `su`, or `supg`) 316 - `none` 317 - 318 319* - `-CtauS` 320 - Scale coefficient for stabilization tau (nondimensional) 321 - `0` 322 - 323 324* - `-wind_type` 325 - Wind type in Advection (`rotation` or `translation`) 326 - `rotation` 327 - 328 329* - `-wind_translation` 330 - Constant wind vector when `-wind_type translation` 331 - `1,0,0` 332 - 333 334* - `-E_wind` 335 - Total energy of inflow wind when `-wind_type translation` 336 - `1E6` 337 - `J` 338 339* - `-bubble_type` 340 - `sphere` (3D) or `cylinder` (2D) 341 - `shpere` 342 - 343 344* - `-bubble_continuity` 345 - `smooth`, `back_sharp`, or `thick` 346 - `smooth` 347 - 348::: 349 350An example of the `rotation` mode can be run with: 351 352``` 353./navierstokes -problem advection -dm_plex_box_faces 10,10,10 -dm_plex_dim 3 -dm_plex_box_lower 0,0,0 -dm_plex_box_upper 8000,8000,8000 -bc_wall 1,2,3,4,5,6 -wall_comps 4 -wind_type rotation -implicit -stab su 354``` 355 356and the `translation` mode with: 357 358``` 359./navierstokes -problem advection -dm_plex_box_faces 10,10,10 -dm_plex_dim 3 -dm_plex_box_lower 0,0,0 -dm_plex_box_upper 8000,8000,8000 -wind_type translation -wind_translation .5,-1,0 -bc_inflow 1,2,3,4,5,6 360``` 361 362### Inviscid Ideal Gas 363 364#### Isentropic Euler vortex 365 366For the Isentropic Vortex problem, the following additional command-line options are available: 367 368:::{list-table} Isentropic Vortex Runtime Options 369:header-rows: 1 370 371* - Option 372 - Description 373 - Default value 374 - Unit 375 376* - `-center` 377 - Location of vortex center 378 - `(lx,ly,lz)/2` 379 - `(m,m,m)` 380 381* - `-units_meter` 382 - 1 meter in scaled length units 383 - `1E-2` 384 - 385 386* - `-units_second` 387 - 1 second in scaled time units 388 - `1E-2` 389 - 390 391* - `-mean_velocity` 392 - Background velocity vector 393 - `(1,1,0)` 394 - 395 396* - `-vortex_strength` 397 - Strength of vortex < 10 398 - `5` 399 - 400 401* - `-c_tau` 402 - Stabilization constant 403 - `0.5` 404 - 405::: 406 407This problem can be run with: 408 409``` 410./navierstokes -problem euler_vortex -dm_plex_box_faces 20,20,1 -dm_plex_box_lower 0,0,0 -dm_plex_box_upper 1000,1000,50 -dm_plex_dim 3 -bc_inflow 4,6 -bc_outflow 3,5 -bc_slip_z 1,2 -mean_velocity .5,-.8,0. 411``` 412 413#### Sod shock tube 414 415For the Shock Tube problem, the following additional command-line options are available: 416 417:::{list-table} Shock Tube Runtime Options 418:header-rows: 1 419 420* - Option 421 - Description 422 - Default value 423 - Unit 424 425* - `-units_meter` 426 - 1 meter in scaled length units 427 - `1E-2` 428 - 429 430* - `-units_second` 431 - 1 second in scaled time units 432 - `1E-2` 433 - 434 435* - `-yzb` 436 - Use YZB discontinuity capturing 437 - `none` 438 - 439 440* - `-stab` 441 - Stabilization method (`none`, `su`, or `supg`) 442 - `none` 443 - 444::: 445 446This problem can be run with: 447 448``` 449./navierstokes -problem shocktube -yzb -stab su -bc_slip_z 3,4 -bc_slip_y 1,2 -bc_wall 5,6 -dm_plex_dim 3 -dm_plex_box_lower 0,0,0 -dm_plex_box_upper 1000,100,100 -dm_plex_box_faces 200,1,1 -units_second 0.1 450``` 451 452### Newtonian viscosity, Ideal Gas 453 454For the Density Current, Channel, and Blasius problems, the following common command-line options are available: 455 456:::{list-table} Newtonian Ideal Gas problems Runtime Options 457:header-rows: 1 458 459* - Option 460 - Description 461 - Default value 462 - Unit 463 464* - `-units_meter` 465 - 1 meter in scaled length units 466 - `1` 467 - 468 469* - `-units_second` 470 - 1 second in scaled time units 471 - `1` 472 - 473 474* - `-units_kilogram` 475 - 1 kilogram in scaled mass units 476 - `1` 477 - 478 479* - `-units_Kelvin` 480 - 1 Kelvin in scaled temperature units 481 - `1` 482 - 483 484* - `-stab` 485 - Stabilization method (`none`, `su`, or `supg`) 486 - `none` 487 - 488 489* - `-c_tau` 490 - Stabilization constant, $c_\tau$ 491 - `0.5` 492 - 493 494* - `-Ctau_t` 495 - Stabilization time constant, $C_t$ 496 - `1.0` 497 - 498 499* - `-Ctau_v` 500 - Stabilization viscous constant, $C_v$ 501 - `36, 60, 128 for degree = 1, 2, 3` 502 - 503 504* - `-Ctau_C` 505 - Stabilization continuity constant, $C_c$ 506 - `1.0` 507 - 508 509* - `-Ctau_M` 510 - Stabilization momentum constant, $C_m$ 511 - `1.0` 512 - 513 514* - `-Ctau_E` 515 - Stabilization energy constant, $C_E$ 516 - `1.0` 517 - 518 519* - `-cv` 520 - Heat capacity at constant volume 521 - `717` 522 - `J/(kg K)` 523 524* - `-cp` 525 - Heat capacity at constant pressure 526 - `1004` 527 - `J/(kg K)` 528 529* - `-g` 530 - Gravitational acceleration 531 - `9.81` 532 - `m/s^2` 533 534* - `-lambda` 535 - Stokes hypothesis second viscosity coefficient 536 - `-2/3` 537 - 538 539* - `-mu` 540 - Shear dynamic viscosity coefficient 541 - `75` 542 - `Pa s` 543 544* - `-k` 545 - Thermal conductivity 546 - `0.02638` 547 - `W/(m K)` 548 549* - `-newtonian_unit_tests` 550 - Developer option to test properties 551 - `false` 552 - boolean 553 554* - `-state_var` 555 - State variables to solve solution with. `conservative` ($\rho, \rho \bm{u}, \rho e$) or `primitive` ($P, \bm{u}, T$) 556 - `conservative` 557 - string 558::: 559 560#### Newtonian Wave 561 562The newtonian wave problem has the following command-line options in addition to the Newtonian Ideal Gas options: 563 564:::{list-table} Newtonian Wave Runtime Options 565:header-rows: 1 566 567* - Option 568 - Description 569 - Default value 570 - Unit 571 572* - `-freestream_riemann` 573 - Riemann solver for boundaries (HLL or HLLC) 574 - `hllc` 575 - 576 577* - `-freestream_velocity` 578 - Freestream velocity vector 579 - `0,0,0` 580 - `m/s` 581 582* - `-freestream_temperature` 583 - Freestream temperature 584 - `288` 585 - `K` 586 587* - `-freestream_pressure` 588 - Freestream pressure 589 - `1.01e5` 590 - `Pa` 591 592* - `-epicenter` 593 - Coordinates of center of perturbation 594 - `0,0,0` 595 - `m` 596 597* - `-amplitude` 598 - Amplitude of the perturbation 599 - `0.1` 600 - 601 602* - `-width` 603 - Width parameter of the perturbation 604 - `0.002` 605 - `m` 606 607::: 608 609This problem can be run with the `newtonianwave.yaml` file via: 610 611``` 612./navierstokes -options_file newtonianwave.yaml 613``` 614 615```{literalinclude} ../../../../../examples/fluids/newtonianwave.yaml 616:language: yaml 617``` 618 619#### Vortex Shedding - Flow past Cylinder 620 621The vortex shedding, flow past cylinder problem has the following command-line options in addition to the Newtonian Ideal Gas options: 622 623:::{list-table} Vortex Shedding Runtime Options 624:header-rows: 1 625 626* - Option 627 - Description 628 - Default value 629 - Unit 630 631* - `-freestream_velocity` 632 - Freestream velocity vector 633 - `0,0,0` 634 - `m/s` 635 636* - `-freestream_temperature` 637 - Freestream temperature 638 - `288` 639 - `K` 640 641* - `-freestream_pressure` 642 - Freestream pressure 643 - `1.01e5` 644 - `Pa` 645 646::: 647 648The initial condition is taken from `-reference_temperature` and `-reference_pressure`. 649To run this problem, first generate a mesh: 650 651```console 652$ make -C examples/fluids/meshes 653``` 654 655Then run by building the executable and running: 656 657```console 658$ make build/fluids-navierstokes 659$ mpiexec -n 6 build/fluids-navierstokes -options_file vortexshedding.yaml 660``` 661 662The vortex shedding period is roughly 6 and this problem runs until time 100 (2000 time steps). 663 664```{literalinclude} ../../../../../examples/fluids/vortexshedding.yaml 665:language: yaml 666``` 667 668#### Density current 669 670The Density Current problem has the following command-line options in addition to the Newtonian Ideal Gas options: 671 672:::{list-table} Density Current Runtime Options 673:header-rows: 1 674 675* - Option 676 - Description 677 - Default value 678 - Unit 679 680* - `-center` 681 - Location of bubble center 682 - `(lx,ly,lz)/2` 683 - `(m,m,m)` 684 685* - `-dc_axis` 686 - Axis of density current cylindrical anomaly, or `(0,0,0)` for spherically symmetric 687 - `(0,0,0)` 688 - 689 690* - `-rc` 691 - Characteristic radius of thermal bubble 692 - `1000` 693 - `m` 694 695* - `-theta0` 696 - Reference potential temperature 697 - `300` 698 - `K` 699 700* - `-thetaC` 701 - Perturbation of potential temperature 702 - `-15` 703 - `K` 704 705* - `-P0` 706 - Atmospheric pressure 707 - `1E5` 708 - `Pa` 709 710* - `-N` 711 - Brunt-Vaisala frequency 712 - `0.01` 713 - `1/s` 714::: 715 716This problem can be run with: 717 718``` 719./navierstokes -problem density_current -dm_plex_box_faces 16,1,8 -degree 1 -dm_plex_box_lower 0,0,0 -dm_plex_box_upper 2000,125,1000 -dm_plex_dim 3 -rc 400. -bc_wall 1,2,5,6 -wall_comps 1,2,3 -bc_slip_y 3,4 -mu 75 720``` 721 722#### Channel flow 723 724The Channel problem has the following command-line options in addition to the Newtonian Ideal Gas options: 725 726:::{list-table} Channel Runtime Options 727:header-rows: 1 728 729* - Option 730 - Description 731 - Default value 732 - Unit 733 734* - `-umax` 735 - Maximum/centerline velocity of the flow 736 - `10` 737 - `m/s` 738 739* - `-theta0` 740 - Reference potential temperature 741 - `300` 742 - `K` 743 744* - `-P0` 745 - Atmospheric pressure 746 - `1E5` 747 - `Pa` 748 749* - `-body_force_scale` 750 - Multiplier for body force (`-1` for flow reversal) 751 - 1 752 - 753::: 754 755This problem can be run with the `channel.yaml` file via: 756 757``` 758./navierstokes -options_file channel.yaml 759``` 760```{literalinclude} ../../../../../examples/fluids/channel.yaml 761:language: yaml 762``` 763 764#### Blasius boundary layer 765 766The Blasius problem has the following command-line options in addition to the Newtonian Ideal Gas options: 767 768:::{list-table} Blasius Runtime Options 769:header-rows: 1 770 771* - Option 772 - Description 773 - Default value 774 - Unit 775 776* - `-velocity_infinity` 777 - Freestream velocity 778 - `40` 779 - `m/s` 780 781* - `-temperature_infinity` 782 - Freestream temperature 783 - `288` 784 - `K` 785 786* - `-temperature_wall` 787 - Wall temperature 788 - `288` 789 - `K` 790 791* - `-delta0` 792 - Boundary layer height at the inflow 793 - `4.2e-3` 794 - `m` 795 796* - `-P0` 797 - Atmospheric pressure 798 - `1.01E5` 799 - `Pa` 800 801* - `-platemesh_refine_height` 802 - Height at which `-platemesh_Ndelta` number of elements should refined into 803 - `5.9E-4` 804 - `m` 805 806* - `-platemesh_Ndelta` 807 - Number of elements to keep below `-platemesh_refine_height` 808 - `45` 809 - 810 811* - `-platemesh_growth` 812 - Growth rate of the elements in the refinement region 813 - `1.08` 814 - 815 816* - `-platemesh_top_angle` 817 - Downward angle of the top face of the domain. This face serves as an outlet. 818 - `5` 819 - `degrees` 820 821* - `-stg_use` 822 - Whether to use stg for the inflow conditions 823 - `false` 824 - 825 826* - `-platemesh_y_node_locs_path` 827 - Path to file with y node locations. If empty, will use mesh warping instead. 828 - `""` 829 - 830 831* - `-n_chebyshev` 832 - Number of Chebyshev terms 833 - `20` 834 - 835 836* - `-chebyshev_` 837 - Prefix for Chebyshev snes solve 838 - 839 - 840 841::: 842 843This problem can be run with the `blasius.yaml` file via: 844 845``` 846./navierstokes -options_file blasius.yaml 847``` 848 849```{literalinclude} ../../../../../examples/fluids/blasius.yaml 850:language: yaml 851``` 852 853#### STG Inflow for Flat Plate 854 855Using the STG Inflow for the blasius problem adds the following command-line options: 856 857:::{list-table} Blasius Runtime Options 858:header-rows: 1 859 860* - Option 861 - Description 862 - Default value 863 - Unit 864 865* - `-stg_inflow_path` 866 - Path to the STGInflow file 867 - `./STGInflow.dat` 868 - 869 870* - `-stg_rand_path` 871 - Path to the STGRand file 872 - `./STGRand.dat` 873 - 874 875* - `-stg_alpha` 876 - Growth rate of the wavemodes 877 - `1.01` 878 - 879 880* - `-stg_u0` 881 - Convective velocity, $U_0$ 882 - `0.0` 883 - `m/s` 884 885* - `-stg_mean_only` 886 - Only impose the mean velocity (no fluctutations) 887 - `false` 888 - 889 890* - `-stg_strong` 891 - Strongly enforce the STG inflow boundary condition 892 - `false` 893 - 894 895* - `-stg_fluctuating_IC` 896 - "Extrude" the fluctuations through the domain as an initial condition 897 - `false` 898 - 899 900::: 901 902This problem can be run with the `blasius.yaml` file via: 903 904``` 905./navierstokes -options_file blasius.yaml -stg_use true 906``` 907 908Note the added `-stg_use true` flag 909This overrides the `stg: use: false` setting in the `blasius.yaml` file, enabling the use of the STG inflow. 910