G12-Sundry-聚焦波造波代码，波浪参数转化

1.Fluent中的UDF造波代码
 

// freakBAN75.c
/*
inlet moving mesh : paddles
outlet damping
fluent udf
*/
#include<stdio.h>
#include"udf.h"
#define N 75
#define Xp 10.0
#define Tp 15.0
#define pi  M_PI
#define g   9.81
#define h   3
#define L0  20
#define L1  30
#define dw  (4.83459-1.17196)/N
#define LW  5
DEFINE_CG_MOTION(ban_moving, dt, cg_vel, cg_omega, time, dtime)
{
	real Sn[75] = { 3.98E-05,7.69E-05,0.000128005,0.000189661,0.000256353,0.000322136,0.0003819,0.000432067,0.000470719,0.000497378,0.000512633,0.000517739,0.000514288,0.000503961,0.000488366,0.000468946,0.000446936,0.000423355,0.000399016,0.000374549,0.000350429,0.000327,0.000304501,0.000283088,0.000262853,0.000243841,0.000226057,0.000209483,0.00019408,0.0001798,0.000166583,0.000154367,0.000143089,0.000132685,0.000123091,0.000114247,0.000106097,9.86E-05,9.17E-05,8.53E-05,7.94E-05,7.40E-05,6.90E-05,6.43E-05,6.00E-05,5.61E-05,5.24E-05,4.91E-05,4.59E-05,4.30E-05,4.03E-05,3.78E-05,3.55E-05,3.34E-05,3.14E-05,2.95E-05,2.78E-05,2.62E-05,2.46E-05,2.32E-05,2.19E-05,2.07E-05,1.96E-05,1.85E-05,1.75E-05,1.66E-05,1.57E-05,1.48E-05,1.41E-05,1.33E-05,1.27E-05,1.20E-05,1.14E-05,1.09E-05,1.03E-05 };
	real k[75] = { 0.232371455840487,0.249786814149107,0.254149656827340,0.274659964303511,0.285474038790606,0.297912164299048,0.313923158777687,0.321324277725442,0.332577950485603,0.345229949359504,0.360209632718581,0.373869700821059,0.386704443897790,0.400905988782487,0.427237180995047,0.431434306971425,0.450262412027362,0.476991100267958,0.482080568587278,0.509089000706641,0.523833411249025,0.542565682494764,0.557217360473855,0.575792813470086,0.608084429196712,0.617619269403596,0.656201811860559,0.679729793272655,0.686334007551193,0.718583444405394,0.748198795453447,0.768706480317268,0.785281925083144,0.830053292544376,0.837797839083065,0.872058923626267,0.908616475332218,0.934559436401070,0.978228744079897,0.998367031019269,1.03032581517445,1.05822982277296,1.08519853842522,1.12335645542617,1.16763283687379,1.17385947072239,1.21893049141042,1.24385049120768,1.29256614364529,1.32825277957608,1.35830589618821,1.39094138295851,1.44126503331233,1.48426098449304,1.52488396596574,1.57111612906692,1.60899440375169,1.63502351616304,1.67799628996486,1.70742018650009,1.74643523019032,1.79116943440899,1.83888251717798,1.89037051470502,1.91890224124605,1.96465747202773,2.01779476405335,2.06291255223232,2.11822833687560,2.13575483404363,2.19861129869732,2.26221147909793,2.29036169441906,2.35191263607766,2.38270458677882 };
	real w, Tr, A, veltemp, T1;
	real f = 0;
	int i;
	for (i = 0; i < N; i++)
	{
		w = pow((g * k[i] * tanh(k[i] * h)), 0.5);
		A = pow((2 * Sn[i] * dw), 0.5);
		Tr = 4 * sinh(k[i] * h) * sinh(k[i] * h) / (sinh(2 * k[i] * h) + (2 * k[i] * h));
		veltemp = w * A / Tr * cos(k[i] * (-Xp) - w * (time - Tp));
		f = f + veltemp;
	}
	w = pow((g * k[0] * tanh(k[0] * h)), 0.5);
	T1 = 2 * pi / w;
	if (time <= (2 * T1))
	{
		cg_vel[0] = f * time / (2 * T1);
	}
	else
	{
		cg_vel[0] = f;
	}
}
/**********************************/
DEFINE_SOURCE(xom_source,c,t,dS,eqn)
{
	real x[ND_ND];
	real source;
	real mu;
	C_CENTROID(x,c,t);
	if(C_U(c,t)<=0)
	{
		mu=2*pow(LW,0.5)*pow(((x[0]-L0)/(L1-L0)),2)*C_R(c,t);
		source=-mu*C_U(c,t);
		dS[eqn]=-mu;
	}
	else
	{
		mu=0;
		source=0;
		dS[eqn]=0;
	}
	return source;
}
/**********************************/
DEFINE_SOURCE(zom_source,c,t,dS,eqn)
{
	real x[ND_ND];
	real source;
	real mu;
	C_CENTROID(x,c,t);
	mu=10*pow(LW,0.5)*pow(((x[0]-L0)/(L1-L0)),2)*C_R(c,t);
	source=-mu*C_W(c,t);
	dS[eqn]=-mu;
	return source;
}
/**********************************/

2.利用Python实现波浪参数转化，并制作GIF动图
focus_wave.py
 

import numpy as np
import math
import matplotlib.pyplot as plt
import matplotlib.animation as animation


class mesh:
    __x_min__ = None
    __x_max__ = None
    __y_min__ = None
    __y_max__ = None
    __x_num__ = None

    def __init__(self, dictionary_arg):
        self.__x_min__ = dictionary_arg['x_min']
        self.__x_max__ = dictionary_arg['x_max']
        self.__y_min__ = dictionary_arg['y_min']
        self.__y_max__ = dictionary_arg['y_max']
        self.__x_num__ = dictionary_arg['x_num']

    def __get_xp__(self):
        num = self.__x_num__ + 1
        xp = np.linspace(self.__x_min__, self.__x_max__, num)
        return xp

    def get_xx(self):
        xp = self.__get_xp__()
        xc = []
        for ii_in in range(len(xp) - 1):
            xm = (xp[ii_in] + xp[ii_in + 1]) / 2.0
            xc.append(xm)
        return xc


class wave_para_converter:
    __water_depth__ = None
    __wave_dw__ = None
    __x_focus__ = None
    __t_focus__ = None
    __wave_sn__ = None
    __wave_number__ = None

    def __init__(self, dictionary_arg):
        self.__water_depth__ = dictionary_arg['water_depth']
        self.__wave_dw__ = dictionary_arg['wave_dw']
        self.__x_focus__ = dictionary_arg['x_focus']
        self.__t_focus__ = dictionary_arg['t_focus']
        self.__wave_sn__ = dictionary_arg['wave_sn']
        self.__wave_number__ = dictionary_arg['wave_number']

    def wave_height(self):
        sn = self.__wave_sn__
        dw = self.__wave_dw__
        wave_height = 2.0 * pow((2.0 * sn * dw), 0.5)
        return wave_height

    def wave_period(self):
        g = 9.81
        k = self.__wave_number__
        h = self.__water_depth__
        wave_omega = pow((g * k * math.tanh(k * h)), 0.5)
        wave_period = (2.0 * math.pi) / wave_omega
        return wave_period

    def wave_phase(self):
        k = self.__wave_number__
        xp = self.__x_focus__
        tp = self.__t_focus__
        omega = (2.0 * math.pi) / (self.wave_period())
        phase = - k * xp + omega * tp
        return phase


class regular_wave:
    __water_depth__ = None
    __wave_height__ = None
    __wave_phase__ = None
    __wave_period__ = None
    __wave_length__ = None

    def __init__(self, dictionary_arg):
        self.__water_depth__ = dictionary_arg['water_depth']
        self.__wave_height__ = dictionary_arg['wave_height']
        self.__wave_phase__ = dictionary_arg['wave_phase']
        self.__wave_period__ = dictionary_arg['wave_period']
        self.__wave_length__ = self.wave_length()

    def wave_speed(self):
        g_mag = 9.81
        wave_number = self.wave_number()
        water_depth = self.__water_depth__
        wave_speed = math.sqrt((g_mag / wave_number) * math.tanh(wave_number * water_depth))
        return wave_speed

    def wave_length(self):
        water_depth = self.__water_depth__
        wave_period = self.__wave_period__
        g = 9.81
        wave_length0 = abs(g) * wave_period * wave_period / (2.0 * math.pi)
        wave_length = wave_length0
        # stokes I
        for ii_in in range(1, 101):
            wave_length = wave_length0 * math.tanh(2.0 * math.pi * water_depth / wave_length)

        return wave_length

    def wave_number(self):
        return (2.0 * math.pi) / (self.wave_length())

    def wave_omega(self):
        return (2.0 * math.pi) / self.__wave_period__

    def wave_surface(self, xx_arg, tt_arg):
        wave_number = self.wave_number()
        wave_omega = self.wave_omega()
        wave_height = self.__wave_height__
        wave_phase = self.__wave_phase__
        theta00 = wave_number * xx_arg - wave_omega * tt_arg + wave_phase
        wave_surface = (wave_height / 2.0) * math.cos(theta00)
        wave_speed = self.wave_speed()
        distance = wave_speed * tt_arg
        if xx_arg > distance:
            wave_surface = 0.0
        return wave_surface


class focus_wave:
    __wave_sn__ = None
    __wave_number__ = None
    __water_depth__ = None
    __wave_dw__ = None
    __x_focus__ = None
    __t_focus__ = None
    __n_focus__ = None

    def __init__(self, dictionary_arg):
        self.__wave_sn__ = dictionary_arg['wave_sn']
        self.__wave_number__ = dictionary_arg['wave_number']
        self.__water_depth__ = dictionary_arg['water_depth']
        self.__wave_dw__ = dictionary_arg['wave_dw']
        self.__x_focus__ = dictionary_arg['x_focus']
        self.__t_focus__ = dictionary_arg['t_focus']
        self.__n_focus__ = dictionary_arg['n_focus']

    def wave_surface(self, xx_arg0, tt_arg0):
        wave_surface = 0.0
        for ii_in in range(self.__n_focus__):
            wpc_dictionary = {
                'water_depth': self.__water_depth__,
                'wave_dw': self.__wave_dw__,
                'x_focus': self.__x_focus__,
                't_focus': self.__t_focus__,
                'wave_sn': self.__wave_sn__[ii_in],
                'wave_number': self.__wave_number__[ii_in]
            }
            my_wpc = wave_para_converter(dictionary_arg=wpc_dictionary)
            regular_wave_dictionary = {
                'water_depth': 3.0,
                'wave_height': my_wpc.wave_height(),
                'wave_period': my_wpc.wave_period(),
                'wave_phase': my_wpc.wave_phase()
            }
            my_regular_wave = regular_wave(dictionary_arg=regular_wave_dictionary)
            wave_surface_tmp = my_regular_wave.wave_surface(xx_arg=xx_arg0, tt_arg=tt_arg0)
            wave_surface = wave_surface_tmp + wave_surface

        return wave_surface

    def write_paras(self):
        my_out_file_name = './waveDict.txt'
        my_out_file = open(my_out_file_name, 'w')
        my_out_file.write('wavePeriods\n')
        my_out_file.write(str(self.__n_focus__) + '\n')
        my_out_file.write('(\n')
        for ii_in in range(self.__n_focus__):
            wpc_dictionary = {
                'water_depth': self.__water_depth__,
                'wave_dw': self.__wave_dw__,
                'x_focus': self.__x_focus__,
                't_focus': self.__t_focus__,
                'wave_sn': self.__wave_sn__[ii_in],
                'wave_number': self.__wave_number__[ii_in]
            }
            my_wpc = wave_para_converter(dictionary_arg=wpc_dictionary)
            my_out_file.write(str(my_wpc.wave_period()) + '\n')
        my_out_file.write(');\n')
        my_out_file.write('waveHeights\n')
        my_out_file.write(str(self.__n_focus__) + '\n')
        my_out_file.write('(\n')
        for ii_in in range(self.__n_focus__):
            wpc_dictionary = {
                'water_depth': self.__water_depth__,
                'wave_dw': self.__wave_dw__,
                'x_focus': self.__x_focus__,
                't_focus': self.__t_focus__,
                'wave_sn': self.__wave_sn__[ii_in],
                'wave_number': self.__wave_number__[ii_in]
            }
            my_wpc = wave_para_converter(dictionary_arg=wpc_dictionary)
            my_out_file.write(str(my_wpc.wave_height()) + '\n')
        my_out_file.write(');\n')
        my_out_file.write('wavePhases\n')
        my_out_file.write(str(self.__n_focus__) + '\n')
        my_out_file.write('(\n')
        for ii_in in range(self.__n_focus__):
            wpc_dictionary = {
                'water_depth': self.__water_depth__,
                'wave_dw': self.__wave_dw__,
                'x_focus': self.__x_focus__,
                't_focus': self.__t_focus__,
                'wave_sn': self.__wave_sn__[ii_in],
                'wave_number': self.__wave_number__[ii_in]
            }
            my_wpc = wave_para_converter(dictionary_arg=wpc_dictionary)
            wave_phase = my_wpc.wave_phase()  # / 180.0 * math.pi
            my_out_file.write(str(wave_phase) + '\n')
        my_out_file.write(');\n')
        my_out_file.write('waveDirs\n')
        my_out_file.write(str(self.__n_focus__) + '\n')
        my_out_file.write('(\n')
        for ii_in in range(self.__n_focus__):
            my_out_file.write('0.0\n')
        my_out_file.write(');\n')
        my_out_file.flush()


if __name__ == '__main__':
    mesh_dictionary = {
        'x_min': 0,
        'x_max': 30,
        'y_min': -0.5,
        'y_max': 0.5,
        'x_num': 100
    }
    my_mesh = mesh(dictionary_arg=mesh_dictionary)
    xx = my_mesh.get_xx()
    focus_dictionary = {'water_depth': 3.0, 'wave_dw': (4.83459 - 1.17196) / 75.0, 'x_focus': 10.0, 't_focus': 15.0,
                        'n_focus': 75,
                        'wave_sn': [3.98E-05, 7.69E-05, 0.000128005, 0.000189661, 0.000256353, 0.000322136, 0.0003819,
                                    0.000432067, 0.000470719, 0.000497378, 0.000512633, 0.000517739, 0.000514288,
                                    0.000503961, 0.000488366, 0.000468946, 0.000446936, 0.000423355, 0.000399016,
                                    0.000374549, 0.000350429, 0.000327, 0.000304501, 0.000283088, 0.000262853,
                                    0.000243841, 0.000226057, 0.000209483, 0.00019408, 0.0001798, 0.000166583,
                                    0.000154367, 0.000143089, 0.000132685, 0.000123091, 0.000114247, 0.000106097,
                                    9.86E-05, 9.17E-05, 8.53E-05, 7.94E-05, 7.40E-05, 6.90E-05, 6.43E-05, 6.00E-05,
                                    5.61E-05, 5.24E-05, 4.91E-05, 4.59E-05, 4.30E-05, 4.03E-05, 3.78E-05, 3.55E-05,
                                    3.34E-05, 3.14E-05, 2.95E-05, 2.78E-05, 2.62E-05, 2.46E-05, 2.32E-05, 2.19E-05,
                                    2.07E-05, 1.96E-05, 1.85E-05, 1.75E-05, 1.66E-05, 1.57E-05, 1.48E-05, 1.41E-05,
                                    1.33E-05, 1.27E-05, 1.20E-05, 1.14E-05, 1.09E-05, 1.03E-05],
                        'wave_number': [0.232371455840487, 0.249786814149107, 0.254149656827340, 0.274659964303511,
                                        0.285474038790606, 0.297912164299048, 0.313923158777687, 0.321324277725442,
                                        0.332577950485603, 0.345229949359504, 0.360209632718581, 0.373869700821059,
                                        0.386704443897790, 0.400905988782487, 0.427237180995047, 0.431434306971425,
                                        0.450262412027362, 0.476991100267958, 0.482080568587278, 0.509089000706641,
                                        0.523833411249025, 0.542565682494764, 0.557217360473855, 0.575792813470086,
                                        0.608084429196712, 0.617619269403596, 0.656201811860559, 0.679729793272655,
                                        0.686334007551193, 0.718583444405394, 0.748198795453447, 0.768706480317268,
                                        0.785281925083144, 0.830053292544376, 0.837797839083065, 0.872058923626267,
                                        0.908616475332218, 0.934559436401070, 0.978228744079897, 0.998367031019269,
                                        1.03032581517445, 1.05822982277296, 1.08519853842522, 1.12335645542617,
                                        1.16763283687379, 1.17385947072239, 1.21893049141042, 1.24385049120768,
                                        1.29256614364529, 1.32825277957608, 1.35830589618821, 1.39094138295851,
                                        1.44126503331233, 1.48426098449304, 1.52488396596574, 1.57111612906692,
                                        1.60899440375169, 1.63502351616304, 1.67799628996486, 1.70742018650009,
                                        1.74643523019032, 1.79116943440899, 1.83888251717798, 1.89037051470502,
                                        1.91890224124605, 1.96465747202773, 2.01779476405335, 2.06291255223232,
                                        2.11822833687560, 2.13575483404363, 2.19861129869732, 2.26221147909793,
                                        2.29036169441906, 2.35191263607766, 2.38270458677882]}
    my_focus = focus_wave(dictionary_arg=focus_dictionary)
    my_focus.write_paras()

    fig, ax = plt.subplots()

    # tt = 15.0
    # phi = []
    # for ii in range(len(xx)):
    #     phi_tmp = my_focus.wave_surface(xx_arg0=xx[ii], tt_arg0=tt)
    #     phi.append(phi_tmp)
    #
    # ax.set_xlim([0.0, 30.0])
    # ax.set_ylim([-0.3, 0.3])
    # line = ax.plot(xx, phi, '-k')
    # plt.show()

    tt_max = 30.0
    tt_step = 0.3
    tt_min = 0.0
    # for tt in np.arange(tt_min, tt_max, tt_step):
    #     plt.cla()
    #     phi = []
    #     for ii in range(len(xx)):
    #         phi_tmp = my_focus.wave_surface(xx_arg0=xx[ii], tt_arg0=tt)
    #         phi.append(phi_tmp)
    #
    #     ax.set_xlim([0.0, 30.0])
    #     ax.set_ylim([-0.3, 0.3])
    #
    #     line = ax.plot(xx, phi, '-k')
    #     plt.pause(0.001)
    line, = ax.plot([], [])
    ax.set_xlim([0.0, 30.0])
    ax.set_ylim([-0.3, 0.3])
    ims = []
    for tt in np.arange(tt_min, tt_max, tt_step):
        phi = []
        for ii in range(len(xx)):
            phi_tmp = my_focus.wave_surface(xx_arg0=xx[ii], tt_arg0=tt)
            phi.append(phi_tmp)

        line, = ax.plot(xx, phi, '-k')
        ims.append([line])

    ani = animation.ArtistAnimation(fig, ims, interval=10, repeat=False)

    writer = animation.PillowWriter(fps=30)
    ani.save('wave.gif', writer=writer)

参考文献