Raw File
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% This file is part of NS2DDV.                                                      %
%                                                                                   %
% Copyright(C) 2011-2018    C. Calgaro  (caterina.calgaro@math.univ-lille1.fr)      %
%                           E. Creusé   (emmanuel.creuse@math.univ-lille1.fr)       %
%                           T. Goudon   (thierry.goudon@inria.fr)                   %
%                           A. Mouton   (alexandre.mouton@math.univ-lille1.fr)      %
%                                                                                   %
% NS2DDV is free software: you can redistribute it and/or modify it under the terms %
% of the GNU General Public License as published by the Free Software Foundation,   %
% either version 3 of the License, or (at your option) any later version.           %
%                                                                                   %
% NS2DDV is distributed in the hope that it will be useful, but WITHOUT ANY         %
% WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A   %
% PARTICULAR PURPOSE. See the GNU General Public License for more details.          %
%                                                                                   %
% You should have received a copy of the GNU General Public License along with      %
% NS2DDV. If not, see <http://www.gnu.org/licenses/>.                               %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function[PARAMETERS] = complete_manual_setup(OLDPARAMETERS)

	% INPUT
	% OLDPARAMETERS		Raw list of simulation parameters obtained from a setup file
	%
	% OUTPUT
	% PARAMETERS		Enhanced list of simulation parameters

	switch OLDPARAMETERS.CV_STUDY
		case {'NONE'}
			PARAMETERS{1} = OLDPARAMETERS;

		case {'IN_SPACE', 'IN_SPACETIME'}
			switch OLDPARAMETERS.MESH.GENERATION
				case {'NS2DDV'}
					switch OLDPARAMETERS.DOMAIN.GEOMETRY
						case {'RECTANGLE'}
							if (numel(OLDPARAMETERS.MESH.NBSEG_X) == numel(OLDPARAMETERS.MESH.NBSEG_Y))
								for kr=1:numel(OLDPARAMETERS.MESH.NBSEG_X)
									PARAMETERS{kr} = OLDPARAMETERS;
									PARAMETERS{kr}.MESH.NBSEG_X = OLDPARAMETERS.MESH.NBSEG_X{kr};
									PARAMETERS{kr}.MESH.NBSEG_Y = OLDPARAMETERS.MESH.NBSEG_Y{kr};
									PARAMETERS{kr}.MESH.NBSEG_X_LIST = OLDPARAMETERS.MESH.NBSEG_X;
									PARAMETERS{kr}.MESH.NBSEG_Y_LIST = OLDPARAMETERS.MESH.NBSEG_Y;
								end
							else
								error('PARAMETERS.MESH.NBSEG_X and PARAMETERS.MESH.NBSEG_Y must have the same number of elements');
							end

						case {'CIRCLE'}
							for kr=1:numel(OLDPARAMETERS.MESH.NBSEG_C)
								PARAMETERS{kr} = OLDPARAMETERS;
								PARAMETERS{kr}.MESH.NBSEG_C = OLDPARAMETERS.MESH.NBSEG_C{kr};
								PARAMETERS{kr}.MESH.NBSEG_C_LIST = OLDPARAMETERS.MESH.NBSEG_C;
							end

						otherwise
							error('Only RECTANGLE and CIRCLE domain geometry are allowed with NS2DDV mesh generation method');
					end

				case {'PDET'}
					for kr=1:numel(OLDPARAMETERS.MESH.H0)
						PARAMETERS{kr} = OLDPARAMETERS;
						PARAMETERS{kr}.MESH.H0 = OLDPARAMETERS.MESH.H0{kr};
						PARAMETERS{kr}.MESH.H0_LIST = OLDPARAMETERS.MESH.H0;
					end

				otherwise
					error('Only NS2DDV and PDET mesh generation methods are allowed in a convergence/analysis study');
			end

		case {'IN_TIME'}
			for kr=1:numel(OLDPARAMETERS.FE.C_STEP_TIME)
				PARAMETERS{kr} = OLDPARAMETERS;
				PARAMETERS{kr}.FE.C_STEP_TIME = OLDPARAMETERS.FE.C_STEP_TIME{kr};
				PARAMETERS{kr}.FE.C_STEP_TIME_LIST = OLDPARAMETERS.FE.C_STEP_TIME;
			end

		case {'IN_REYNOLDS'}
			for kr=1:numel(OLDPARAMETERS.PHYSICAL.RE)
				PARAMETERS{kr} = OLDPARAMETERS;
				PARAMETERS{kr}.PHYSICAL.RE = OLDPARAMETERS.PHYSICAL.RE{kr};
				PARAMETERS{kr}.PHYSICAL.RE_LIST = OLDPARAMETERS.PHYSICAL.RE;
				PARAMETERS{kr}.FE.ROT_CORRECTION = OLDPARAMETERS.FE.ROT_CORRECTION{kr};
				PARAMETERS{kr}.FE.ROT_CORRECTION_LIST = OLDPARAMETERS.FE.ROT_CORRECTION;
			end

		otherwise
			error('Wrong value for CV_STUDY');
	end


	for kr=1:numel(PARAMETERS)
		[~,PARAMETERS{kr}.HOSTNAME] = system('hostname');

	    switch PARAMETERS{kr}.MODEL
    	    case {'NSDV'}
    	        % Debug mode
    	        % TRANSPORT_ONLY : only solve transport with FV scheme, the velocity is computed with the exact solution
    	        % STOKES_ONLY : only solve Stokes with FE scheme, the density is computed with the exact solution
    	        % Other values : Navier-Stokes resolution with the full hybrid method
    	        PARAMETERS{kr}.DEBUG_MODE = '';
    	    case {'NS'}
    	        PARAMETERS{kr}.DEBUG_MODE = '';
    	    otherwise
				error(sprintf('Unknown model\n'));
    	end

		PARAMETERS{kr}.TESTCASE_ANALYTIC = (max(strcmp(PARAMETERS{kr}.TESTCASE, {'EXAC', 'EXACNEU'})) == 1);

		switch PARAMETERS{kr}.TESTCASE
			case {'CAEN', 'RTIN', 'DROP', 'EXAC'}
				PARAMETERS{kr}.FE.FIXPRESSURE = true;
				PARAMETERS{kr}.FE.BC_NATURAL = false;
				PARAMETERS{kr}.FE.BC_STRESS_BALANCE = false;
			case {'EXACNEU'}
				if (strcmp(PARAMETERS{kr}.FE.BC_LEFT_UX, 'NEUMANN') && ...
									strcmp(PARAMETERS{kr}.FE.BC_LEFT_UY, 'NEUMANN'))
					PARAMETERS{kr}.FE.FIXPRESSURE = false;
					PARAMETERS{kr}.FE.BC_NATURAL = false;
					PARAMETERS{kr}.FE.BC_STRESS_BALANCE = false;
				else
					error('Problem with the values PARAMETERS.FE.BC_LEFT_U*');
				end
			case {'POIS', 'NONA', 'GTPSI'}
				if (strcmp(PARAMETERS{kr}.FE.BC_RIGHT_UX, 'NEUMANN') && ...
									strcmp(PARAMETERS{kr}.FE.BC_RIGHT_UY, 'NEUMANN'))
					PARAMETERS{kr}.FE.FIXPRESSURE = false;
					PARAMETERS{kr}.FE.BC_NATURAL = false;
					PARAMETERS{kr}.FE.BC_STRESS_BALANCE = false;
				elseif (strcmp(PARAMETERS{kr}.FE.BC_RIGHT_UX, 'NATURAL') && ...
									strcmp(PARAMETERS{kr}.FE.BC_RIGHT_UY, 'NATURAL'))
					PARAMETERS{kr}.FE.FIXPRESSURE = true;
					PARAMETERS{kr}.FE.BC_NATURAL = true;
					PARAMETERS{kr}.FE.BC_STRESS_BALANCE = false;
				elseif (strcmp(PARAMETERS{kr}.FE.BC_RIGHT_UX, 'STRESS_BALANCE') && ...
									strcmp(PARAMETERS{kr}.FE.BC_RIGHT_UY, 'STRESS_BALANCE'))
					PARAMETERS{kr}.FE.FIXPRESSURE = false;
					PARAMETERS{kr}.FE.BC_NATURAL = false;
					PARAMETERS{kr}.FE.BC_STRESS_BALANCE = true;
				else
					error('Problem with the values PARAMETERS.FE.BC_RIGHT_U*');
				end
				if max(strcmp(PARAMETERS{kr}.TESTCASE, {'POIS', 'GTPSI'}))
					PARAMETERS{kr}.PHYSICAL.CHAR_VELOCITY = PARAMETERS{kr}.PHYSICAL.INLET_VELOCITY;
				elseif strcmp(PARAMETERS{kr}.TESTCASE, 'NONA')
					switch PARAMETERS{kr}.DOMAIN.GEOMETRY
						case {'RECTANGLE'}
							yinlet = PARAMETERS{kr}.DOMAIN.YMIN;
						case {'STEP'}
							if strcmp(PARAMETERS{kr}.DOMAIN.STEP_POSITION, 'LEFT')
								yinlet = PARAMETERS{kr}.DOMAIN.YSTEP;
							elseif strcmp(PARAMETERS{kr}.DOMAIN.STEP_POSITION, 'RIGHT')
								yinlet = PARAMETERS{kr}.DOMAIN.YMIN;
							else
								error('Wrong value of PARAMETERS.DOMAIN.STEP_POSITION');
							end
						case {'DIHEDRON'}
							yinlet = PARAMETERS{kr}.DOMAIN.YAXIS_ANGLE1;
						otherwise
							error('Unknown domain geometry');
					end
					eta = (PARAMETERS{kr}.DOMAIN.YMAX-yinlet) * sqrt(PARAMETERS{kr}.PHYSICAL.RE / ...
								(PARAMETERS{kr}.DOMAIN.XMIN-PARAMETERS{kr}.PHYSICAL.XSTART));
					PARAMETERS{kr}.PHYSICAL.CHAR_VELOCITY = ...
						sqrt((PARAMETERS{kr}.DOMAIN.XMIN-PARAMETERS{kr}.PHYSICAL.XSTART) / ...
						PARAMETERS{kr}.PHYSICAL.RE)/(PARAMETERS{kr}.DOMAIN.YMAX-yinlet)*fBlasius(eta);
				end

			otherwise
				error('Unknown test case');
		end

		if max(strcmp(PARAMETERS{kr}.TESTCASE, {'NONA'}))
			% The initialization of the test case is not divergence-free
			% We must count the number of time iterations on which the FE scheme must ensure the divergence-free
			% condition (ex: replace 'BDF2_PROJ' by 'BDF2')
			switch PARAMETERS{kr}.FE.SCHEME
				case{'BDF2'}
					% The divergence-free constrain is a part of the scheme
					PARAMETERS{kr}.PHYSICAL.FORCE_INCOMP = 0;
				case {'BDF2_PROJ', 'BDF2_PROJ_VAR', 'BDF2_PROJ_VAR_A'}
					% The divergence-free constrain is not insured
					% --> Replace by 'BDF2' on the 2 first time iterations
					PARAMETERS{kr}.PHYSICAL.FORCE_INCOMP = 2;
				otherwise
					error('Unknown FE scheme');
			end
		end

		% Complete the physical parameters
		switch PARAMETERS{kr}.TESTCASE
			case {'RTIN'}
				% Provide the Atwood number
				%%%%%%%%%%%%%%%%%%%%%%%%%%%
				PARAMETERS{kr}.PHYSICAL.ATWOOD = ...
						(PARAMETERS{kr}.PHYSICAL.RHOMAX-PARAMETERS{kr}.PHYSICAL.RHOMIN) / ...
						(PARAMETERS{kr}.PHYSICAL.RHOMAX+PARAMETERS{kr}.PHYSICAL.RHOMIN);

				% Provide the characteristic size of the domain
				%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
				PARAMETERS{kr}.PHYSICAL.RHOD = PARAMETERS{kr}.DOMAIN.XMAX-PARAMETERS{kr}.DOMAIN.XMIN;

				% Provide the amplitude of the initial perturbation
				%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
				if (PARAMETERS{kr}.PHYSICAL.ATWOOD < 0.5)
					PARAMETERS{kr}.PHYSICAL.RHOETA = 0.01;
				else
					PARAMETERS{kr}.PHYSICAL.RHOETA = 0.1;
				end

				% Provide the amplitude of the viscosity
				%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
				if strcmp(PARAMETERS{kr}.PHYSICAL.VISCOSITY, 'CONST')
					PARAMETERS{kr}.PHYSICAL.MU = ...
							PARAMETERS{kr}.PHYSICAL.RHOMIN*(PARAMETERS{kr}.PHYSICAL.RHOD^1.5) * ...
							sqrt(PARAMETERS{kr}.PHYSICAL.GRAVITY)/PARAMETERS{kr}.PHYSICAL.RE;
				else
					error(sprintf('Only constant viscosity can be handled for a Rayleigh-Taylor instability\n'));
				end

				% Provide the boundary conditions
				%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
				% 'A' : No periodic boundary conditions on left and right bounds
				% 'P' : Periodic boundary conditions on left and right bounds
				PARAMETERS{kr}.CLROG = 'A';
				PARAMETERS{kr}.CLROD = 'A';

			case {'DROP'}
				PARAMETERS{kr}.PHYSICAL.MU = 1.;

				% Provide the boundary conditions
				%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
				% 'A' : No periodic boundary conditions on left and right bounds
				% 'P' : Periodic boundary conditions on left and right bounds
				PARAMETERS{kr}.CLROG = 'A';
				PARAMETERS{kr}.CLROD = 'A';

			case {'EXAC', 'EXACNEU'}%, 'EXAC2'}
				PARAMETERS{kr}.PHYSICAL.MU = 1.;

				% Provide the boundary conditions
				%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
				% 'A' : No periodic boundary conditions on left and right bounds
				% 'P' : Periodic boundary conditions on left and right bounds
				PARAMETERS{kr}.CLROG = 'A';
				PARAMETERS{kr}.CLROD = 'A';

			case {'CAEN'}
				PARAMETERS{kr}.PHYSICAL.MU = 1.;

				% Provide the boundary conditions
				%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
				% 'A' : No periodic boundary conditions on left and right bounds
				% 'P' : Periodic boundary conditions on left and right bounds
				PARAMETERS{kr}.CLROG = 'A';
				PARAMETERS{kr}.CLROD = 'A';

			case {'NONA'}

			case {'POIS'}

			case {'GTPSI'}

			otherwise
				error(sprintf('Unknown test case.\n'));
		end

		switch PARAMETERS{kr}.MODEL
			case {'NS','NSDV'}
				if (strcmp(PARAMETERS{kr}.MESH.GENERATION, 'NS2DDV') ...
									&& strcmp(PARAMETERS{kr}.TESTCASE, {'NONA'}))
					PARAMETERS{kr}.FE.FIT_STEP_TIME = 'NO';
					switch PARAMETERS{kr}.DOMAIN.GEOMETRY
						case {'RECTANGLE'}
							dx = (PARAMETERS{kr}.DOMAIN.XMAX-PARAMETERS{kr}.DOMAIN.XMIN) / ...
											PARAMETERS{kr}.MESH.NBSEG_X;
						case {'DIHEDRON'}
							dx = 1./PARAMETERS{kr}.MESH.NBSEG_PERUNIT_X;
						otherwise
							error('Problem in complete_manual_setup');
					end

					switch PARAMETERS{kr}.DOMAIN.GEOMETRY
						case {'RECTANGLE'}
							hl = PARAMETERS{kr}.DOMAIN.YMAX-PARAMETERS{kr}.DOMAIN.YMIN;
							hr = PARAMETERS{kr}.DOMAIN.YMAX-PARAMETERS{kr}.DOMAIN.YMIN;
						case {'STEP'}
							if strcmp(PARAMETERS{kr}.DOMAIN.STEP_POSITION, 'LEFT')
								hl = PARAMETERS{kr}.DOMAIN.YMAX-PARAMETERS{kr}.DOMAIN.YSTEP;
								hr = PARAMETERS{kr}.DOMAIN.YMAX-PARAMETERS{kr}.DOMAIN.YMIN;
							elseif strcmp(PARAMETERS{kr}.DOMAIN.STEP_POSITION, 'RIGHT')
								yinlet = PARAMETERS{kr}.DOMAIN.YMIN;
								hl = PARAMETERS{kr}.DOMAIN.YMAX-PARAMETERS{kr}.DOMAIN.YMIN;
								hr = PARAMETERS{kr}.DOMAIN.YMAX-PARAMETERS{kr}.DOMAIN.YSTEP;
							else
								error('Wrong value of PARAMETERS.DOMAIN.STEP_POSITION');
							end
						case {'DIHEDRON'}
							hl = PARAMETERS{kr}.DOMAIN.YMAX-PARAMETERS{kr}.DOMAIN.YAXIS_ANGLE1;
							hr = PARAMETERS{kr}.DOMAIN.YMAX-PARAMETERS{kr}.DOMAIN.YMIN;
						otherwise
							error('Unknown domain geometry');
					end
					PARAMETERS{kr}.FE.CMAX_STEP_TIME = dx/(PARAMETERS{kr}.PHYSICAL.CHAR_VELOCITY*(hl/hr));
					PARAMETERS{kr}.FE.CMIN_STEP_TIME = 0.;
					PARAMETERS{kr}.FE.ALPHAMAX_STEP_TIME = 0.;
					PARAMETERS{kr}.FE.ALPHAMIN_STEP_TIME = 0.;
					PARAMETERS{kr}.FE.STEP_TIME_MODIFIED = true;
				else
					PARAMETERS{kr}.FE.STEP_TIME_MODIFIED = false;
				end
			otherwise
				error('Unknown model');
		end

		% Create and/or test the output directory
		%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
		if ~strcmp(PARAMETERS{kr}.CV_STUDY, 'NONE')
			% Provide an output directory per run
			PARAMETERS{kr}.OUTPUT.DIRECTORY_NAME = strcat(OLDPARAMETERS.OUTPUT.DIRECTORY_NAME, '/', ...
						OLDPARAMETERS.MODEL, '_', OLDPARAMETERS.TESTCASE, '_', sprintf('%d', kr));
		end
		PARAMETERS{kr}.OUTPUT.DIRECTORY_NAME = newdir(PARAMETERS{kr}.OUTPUT.DIRECTORY_NAME, './RESULTS');
		expr = strcat(PARAMETERS{kr}.OUTPUT.DIRECTORY_NAME, '/', PARAMETERS{kr}.OUTPUT.FILE_NAME, '_*');
		if ~isfield(PARAMETERS{kr}, 'RESTART')
			% The current run is not initialized with restart_ns2ddv
			testdir(expr);
		else
			if ~PARAMETERS{kr}.RESTART
				% The current run is not initialized with restart_ns2ddv
				testdir(expr);
			end
		end
	
		% Create and/or test the backup directory
		%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
		PARAMETERS{kr}.BACKUP.DIRECTORY_NAME = newdir(PARAMETERS{kr}.BACKUP.DIRECTORY_NAME, './BACKUP');
		expr = strcat(PARAMETERS{kr}.BACKUP.DIRECTORY_NAME, '/backup_*');
		if ~isfield(PARAMETERS{kr}, 'RESTART')
			% The current run is not initialized with restart_ns2ddv
			testdir(expr);
		else
			if ~PARAMETERS{kr}.RESTART
				% The current run is not initialized with restart_ns2ddv
				testdir(expr);
			end
		end

		if ~strcmp(OLDPARAMETERS.OUTPUT.XDISPLAY, 'NONE')
			PARAMETERS{kr}.OUTPUT.XDISPLAY_LIST = unique(sort(OLDPARAMETERS.OUTPUT.XDISPLAY_LIST));
		else
			PARAMETERS{kr}.OUTPUT.XDISPLAY_LIST = {};
		end
		if ~strcmp(OLDPARAMETERS.OUTPUT.XDISPLAY_SAVE, 'NONE')
			PARAMETERS{kr}.OUTPUT.XDISPLAY_SAVE_OPTDIAGS = unique(sort(OLDPARAMETERS.OUTPUT.XDISPLAY_SAVE_OPTDIAGS));
		else
			PARAMETERS{kr}.OUTPUT.XDISPLAY_SAVE_OPTDIAGS = {};
		end
	end

end
back to top