From e313a4775db4a86a11bd9e9456a7cbafd65fbcbe Mon Sep 17 00:00:00 2001 From: James E McClure Date: Mon, 13 Jul 2015 11:17:56 -0400 Subject: [PATCH] Updated component algorithms in TwoPhase.h --- common/TwoPhase.h | 169 ++++++++++++++++++++++++++++++++++++++++++---- 1 file changed, 155 insertions(+), 14 deletions(-) diff --git a/common/TwoPhase.h b/common/TwoPhase.h index 38d87931..0baa1210 100644 --- a/common/TwoPhase.h +++ b/common/TwoPhase.h @@ -5,7 +5,7 @@ #include "analysis/analysis.h" #include -#define BLOB_AVG_COUNT 27 +#define BLOB_AVG_COUNT 28 // Array access for averages defined by the following #define VOL 0 @@ -23,18 +23,19 @@ #define VX 12 #define VY 13 #define VZ 14 -#define VWNX 15 -#define VWNY 16 -#define VWNZ 17 -#define VWNSX 18 -#define VWNSY 19 -#define VWNSZ 20 -#define GWNXX 21 -#define GWNYY 22 -#define GWNZZ 23 -#define GWNXY 24 -#define GWNXZ 25 -#define GWNYZ 26 +#define VSQ 15 +#define VWNX 16 +#define VWNY 17 +#define VWNZ 18 +#define VWNSX 19 +#define VWNSY 20 +#define VWNSZ 21 +#define GWNXX 22 +#define GWNYY 23 +#define GWNZZ 24 +#define GWNXY 25 +#define GWNXZ 26 +#define GWNYZ 27 class TwoPhase{ @@ -74,6 +75,8 @@ class TwoPhase{ // CSV / text file where time history of averages is saved FILE *TIMELOG; + FILE *NWPLOG; + FILE *WPLOG; public: //........................................................................... @@ -231,6 +234,16 @@ public: fprintf(TIMELOG,"trawn trJwn trRwn\n"); // trimmed curvature for wn surface //fprintf(TIMELOG,"--------------------------------------------------------------------------------------\n"); } + + NWPLOG = fopen("components.NWP.tcat","a+"); + fprintf(NWPLOG,"time label vol pn awn ans Jwn Kwn lwns cwns "); + fprintf(NWPLOG,"vx vy vz vwnx vwny vwnz vwnsx vwnsy vwnsz vsq "); + fprintf(NWPLOG,"Gwnxx Gwnyy Gwnzz Gwnxy Gwnxz Gwnyz\n"); + + WPLOG = fopen("components.WP.tcat","a+"); + fprintf(WPLOG,"time label vol pw awn ans Jwn Kwn lwns cwns "); + fprintf(WPLOG,"vx vy vz vwnx vwny vwnz vwnsx vwnsy vwnsz vsq "); + fprintf(WPLOG,"Gwnxx Gwnyy Gwnzz Gwnxy Gwnxz Gwnyz\n"); } } ~TwoPhase(){ @@ -250,7 +263,7 @@ public: void PrintAll(int timestep); int GetCubeLabel(int i, int j, int k, IntArray &BlobLabel); void SortBlobs(); - + void PrintComponents(int timestep); }; void TwoPhase::ColorToSignedDistance(double Beta, DoubleArray &ColorData, DoubleArray &DistData){ @@ -587,6 +600,9 @@ void TwoPhase::ComponentAverages(){ ComponentAverages_NWP(VX,LabelNWP) += 0.125*Vel_x(n); ComponentAverages_NWP(VY,LabelNWP) += 0.125*Vel_y(n); ComponentAverages_NWP(VZ,LabelNWP) += 0.125*Vel_z(n); + // twice the kinetic energy + ComponentAverages_NWP(VSQ,LabelNWP) += 0.125*(Vel_x(n)*Vel_x(n)+Vel_y(n)*Vel_y(n)+Vel_z(n)*Vel_z(n)); + // volume the for pressure averaging excludes the interfacial region if (DelPhi(n) < 1e-4 ){ ComponentAverages_NWP(TRIMVOL,LabelNWP) += 0.125; @@ -600,6 +616,9 @@ void TwoPhase::ComponentAverages(){ ComponentAverages_WP(VX,LabelWP) += 0.125*Vel_x(n); ComponentAverages_WP(VY,LabelWP)+= 0.125*Vel_y(n); ComponentAverages_WP(VZ,LabelWP) += 0.125*Vel_z(n); + // twice the kinetic energy + ComponentAverages_WP(VSQ,LabelWP) += 0.125*(Vel_x(n)*Vel_x(n)+Vel_y(n)*Vel_y(n)+Vel_z(n)*Vel_z(n)); + // volume the for pressure averaging excludes the interfacial region if (DelPhi(n) < 1e-4){ ComponentAverages_WP(TRIMVOL,LabelWP) += 0.125; @@ -710,6 +729,66 @@ void TwoPhase::ComponentAverages(){ } } } + + // Globally reduce the non-wetting phase averages + for (int b=0; b 0.0){ + double Vn,pn,awn,ans,Jwn,Kwn,lwns,cwns,vsq; + + Vn = ComponentAverages_NWP(VOL,b); + awn = ComponentAverages_NWP(AWN,b); + ans = ComponentAverages_NWP(ANS,b); + van(0) = ComponentAverages_NWP(VX,b)/Vn; + van(1) = ComponentAverages_NWP(VY,b)/Vn; + van(2) = ComponentAverages_NWP(VZ,b)/Vn; + vsq = ComponentAverages_NWP(VSQ,b)/Vn; + + if (ComponentAverages_NWP(TRIMVOL,b) > 0.0){ + pn = ComponentAverages_NWP(PRS,b)/ComponentAverages_NWP(TRIMVOL,b); + } + else pn = 0.0; + + if (awn != 0.0){ + Jwn = ComponentAverages_NWP(JWN,b)/awn; + Kwn = ComponentAverages_NWP(KWN,b)/awn; + vawn(0) = ComponentAverages_NWP(VWNSX,b)/awn; + vawn(1) = ComponentAverages_NWP(VWNSY,b)/awn; + vawn(2) = ComponentAverages_NWP(VWNSZ,b)/awn; + } + else Jwn=Kwn=0.0; + + lwns = ComponentAverages_NWP(LWNS,b); + if (lwns != 0.0){ + cwns = ComponentAverages_NWP(CWNS,b)/lwns; + vawns(0) = ComponentAverages_NWP(VWNSX,b)/lwns; + vawns(1) = ComponentAverages_NWP(VWNSY,b)/lwns; + vawns(2) = ComponentAverages_NWP(VWNSZ,b)/lwns; + } + else cwns=0.0; + + ComponentAverages_NWP(PRS,b) = pn; + ComponentAverages_NWP(VNX,b) = van(0); + ComponentAverages_NWP(VNY,b) = van(1); + ComponentAverages_NWP(VNZ,b) = van(2); + ComponentAverages_NWP(VSQ,b) = vsq; + + ComponentAverages_NWP(JWN,b) = Jwn; + ComponentAverages_NWP(KWN,b) = Kwn; + ComponentAverages_NWP(VWNX,b) = vawn(0); + ComponentAverages_NWP(VWNY,b) = vawn(1); + ComponentAverages_NWP(VWNZ,b) = vawn(2); + + ComponentAverages_NWP(CWNS,b) = cwns; + ComponentAverages_NWP(VWNSX,b) = vawns(0); + ComponentAverages_NWP(VWNSY,b) = vawns(1); + ComponentAverages_NWP(VWNSZ,b) = vawns(2); + + } + } } @@ -830,6 +909,68 @@ void TwoPhase::PrintAll(int timestep){ } } +void TwoPhase::PrintComponents(int timestep){ + if (Dm.rank==0){ + for (int b=0; b