Commit 3f9a316d authored by breuera's avatar breuera

Merge branch 'master' of github.com:TUM-I5/SWE

parents 2bd1c3d8 aca161d1
<?xml version="1.0" encoding="utf-8"?>
<Project DefaultTargets="Build" ToolsVersion="4.0" xmlns="http://schemas.microsoft.com/developer/msbuild/2003">
<ItemGroup Label="ProjectConfigurations">
<ProjectConfiguration Include="Debug|Win32">
<Configuration>Debug</Configuration>
<Platform>Win32</Platform>
</ProjectConfiguration>
<ProjectConfiguration Include="Release|Win32">
<Configuration>Release</Configuration>
<Platform>Win32</Platform>
</ProjectConfiguration>
</ItemGroup>
<PropertyGroup Label="Globals">
<ProjectGuid>{88B845B4-569C-49D1-9C07-967D2E58EAC9}</ProjectGuid>
<RootNamespace>SWEOpenGL</RootNamespace>
</PropertyGroup>
<Import Project="$(VCTargetsPath)\Microsoft.Cpp.Default.props" />
<PropertyGroup Condition="'$(Configuration)|$(Platform)'=='Debug|Win32'" Label="Configuration">
<ConfigurationType>Application</ConfigurationType>
<UseDebugLibraries>true</UseDebugLibraries>
<CharacterSet>MultiByte</CharacterSet>
</PropertyGroup>
<PropertyGroup Condition="'$(Configuration)|$(Platform)'=='Release|Win32'" Label="Configuration">
<ConfigurationType>Application</ConfigurationType>
<UseDebugLibraries>false</UseDebugLibraries>
<WholeProgramOptimization>true</WholeProgramOptimization>
<CharacterSet>MultiByte</CharacterSet>
</PropertyGroup>
<Import Project="$(VCTargetsPath)\Microsoft.Cpp.props" />
<ImportGroup Label="ExtensionSettings">
<Import Project="$(VCTargetsPath)\BuildCustomizations\CUDA 4.2.props" />
</ImportGroup>
<ImportGroup Label="PropertySheets" Condition="'$(Configuration)|$(Platform)'=='Debug|Win32'">
<Import Project="$(UserRootDir)\Microsoft.Cpp.$(Platform).user.props" Condition="exists('$(UserRootDir)\Microsoft.Cpp.$(Platform).user.props')" Label="LocalAppDataPlatform" />
</ImportGroup>
<ImportGroup Label="PropertySheets" Condition="'$(Configuration)|$(Platform)'=='Release|Win32'">
<Import Project="$(UserRootDir)\Microsoft.Cpp.$(Platform).user.props" Condition="exists('$(UserRootDir)\Microsoft.Cpp.$(Platform).user.props')" Label="LocalAppDataPlatform" />
</ImportGroup>
<PropertyGroup Label="UserMacros" />
<PropertyGroup Condition="'$(Configuration)|$(Platform)'=='Release|Win32'">
<IncludePath>C:\Program Files\NVIDIA GPU Computing Toolkit\CUDA\v4.2\include;$(IncludePath)</IncludePath>
</PropertyGroup>
<PropertyGroup Condition="'$(Configuration)|$(Platform)'=='Release|Win32'">
<ReferencePath>C:\Program Files\NVIDIA GPU Computing Toolkit\CUDA\v4.2\lib;$(ReferencePath)</ReferencePath>
<LibraryPath>$(LibraryPath)</LibraryPath>
</PropertyGroup>
<ItemDefinitionGroup Condition="'$(Configuration)|$(Platform)'=='Debug|Win32'">
<ClCompile>
<WarningLevel>Level3</WarningLevel>
<Optimization>Disabled</Optimization>
<AdditionalIncludeDirectories>SDL\include;$(CudaToolkitIncludeDir)</AdditionalIncludeDirectories>
<RuntimeLibrary>MultiThreadedDLL</RuntimeLibrary>
</ClCompile>
<Link>
<GenerateDebugInformation>true</GenerateDebugInformation>
<AdditionalLibraryDirectories>SDL\lib;%(AdditionalLibraryDirectories)</AdditionalLibraryDirectories>
<AdditionalDependencies>opengl32.lib;glu32.lib;sdl.lib;sdlmain.lib;cudart.lib;%(AdditionalDependencies)</AdditionalDependencies>
</Link>
<CudaCompile>
<Runtime>MD</Runtime>
<CodeGeneration>compute_20,sm_20</CodeGeneration>
</CudaCompile>
<PostBuildEvent>
<Command>xcopy /y "$(ProjectDir)SDL\lib\SDL.dll" "$(OutDir)" </Command>
</PostBuildEvent>
</ItemDefinitionGroup>
<ItemDefinitionGroup Condition="'$(Configuration)|$(Platform)'=='Release|Win32'">
<ClCompile>
<WarningLevel>Level3</WarningLevel>
<Optimization>MaxSpeed</Optimization>
<FunctionLevelLinking>true</FunctionLevelLinking>
<IntrinsicFunctions>true</IntrinsicFunctions>
<AdditionalIncludeDirectories>SDL\include;$(CudaToolkitIncludeDir)</AdditionalIncludeDirectories>
</ClCompile>
<Link>
<GenerateDebugInformation>true</GenerateDebugInformation>
<EnableCOMDATFolding>true</EnableCOMDATFolding>
<OptimizeReferences>true</OptimizeReferences>
<AdditionalDependencies>opengl32.lib;glu32.lib;sdl.lib;sdlmain.lib;cudart.lib;%(AdditionalDependencies)</AdditionalDependencies>
<AdditionalLibraryDirectories>SDL\lib;%(AdditionalLibraryDirectories)</AdditionalLibraryDirectories>
<SubSystem>Windows</SubSystem>
</Link>
<PostBuildEvent>
<Command>xcopy /y $(ProjectDir)SDL\lib\SDL.dll $(OutDir) </Command>
</PostBuildEvent>
<ProjectReference>
<UseLibraryDependencyInputs>true</UseLibraryDependencyInputs>
</ProjectReference>
</ItemDefinitionGroup>
<ItemGroup>
<ClInclude Include="src\opengl\camera.h" />
<ClInclude Include="src\opengl\controller.h" />
<ClInclude Include="src\opengl\shader.h" />
<ClInclude Include="src\opengl\simulation.h" />
<ClInclude Include="src\opengl\visualization.h" />
<ClInclude Include="src\scenarios\SWE_simple_scenarios_vis.h" />
<ClInclude Include="src\scenarios\SWE_VisInfo.h" />
<ClInclude Include="src\scenarios\SWE_VtkScenario.h" />
<ClInclude Include="src\scenarios\SWE_VtkScenarioVisInfo.h" />
</ItemGroup>
<ItemGroup>
<ClCompile Include="src\examples\swe_opengl.cpp" />
<ClCompile Include="src\opengl\camera.cpp" />
<ClCompile Include="src\opengl\controller.cpp" />
<ClCompile Include="src\opengl\shader.cpp" />
<ClCompile Include="src\opengl\visualization.cpp" />
<ClCompile Include="src\scenarios\SWE_VtkScenario.cpp" />
<ClCompile Include="src\scenarios\SWE_VtkScenarioVisInfo.cpp" />
<ClCompile Include="src\SWE_Block.cpp" />
<ClCompile Include="src\SWE_WavePropagationBlock.cpp" />
</ItemGroup>
<ItemGroup>
<CudaCompile Include="src\opengl\simulation.cu">
<FileType>Document</FileType>
</CudaCompile>
<CudaCompile Include="src\SWE_BlockCUDA.cu" />
<CudaCompile Include="src\SWE_BlockCUDA_kernels.cu" />
<CudaCompile Include="src\SWE_RusanovBlockCUDA.cu" />
<CudaCompile Include="src\SWE_RusanovBlockCUDA_kernels.cu" />
<CudaCompile Include="src\SWE_WavePropagationBlockCuda.cu" />
<CudaCompile Include="src\SWE_WavePropagationBlockCuda_kernels.cu" />
</ItemGroup>
<ItemGroup>
<None Include="src\SWE_Block.hh" />
<None Include="src\SWE_BlockCUDA.hh" />
<None Include="src\SWE_BlockCUDA_kernels.hh" />
<None Include="src\SWE_RusanovBlockCUDA.hh" />
<None Include="src\SWE_RusanovBlockCUDA_kernels.hh" />
<None Include="src\SWE_WavePropagationBlock.hh" />
<None Include="src\SWE_WavePropagationBlockCuda.hh" />
<None Include="src\SWE_WavePropagationBlockCuda_kernels.hh" />
</ItemGroup>
<Import Project="$(VCTargetsPath)\Microsoft.Cpp.targets" />
<ImportGroup Label="ExtensionTargets">
<Import Project="$(VCTargetsPath)\BuildCustomizations\CUDA 4.2.targets" />
</ImportGroup>
</Project>
\ No newline at end of file
<?xml version="1.0" encoding="utf-8"?>
<Project ToolsVersion="4.0" xmlns="http://schemas.microsoft.com/developer/msbuild/2003">
<ItemGroup>
<Filter Include="Quelldateien">
<UniqueIdentifier>{4FC737F1-C7A5-4376-A066-2A32D752A2FF}</UniqueIdentifier>
<Extensions>cpp;c;cc;cxx;def;odl;idl;hpj;bat;asm;asmx</Extensions>
</Filter>
<Filter Include="Headerdateien">
<UniqueIdentifier>{93995380-89BD-4b04-88EB-625FBE52EBFB}</UniqueIdentifier>
<Extensions>h;hpp;hxx;hm;inl;inc;xsd</Extensions>
</Filter>
<Filter Include="Ressourcendateien">
<UniqueIdentifier>{67DA6AB6-F800-4c08-8B7A-83BB121AAD01}</UniqueIdentifier>
<Extensions>rc;ico;cur;bmp;dlg;rc2;rct;bin;rgs;gif;jpg;jpeg;jpe;resx;tiff;tif;png;wav;mfcribbon-ms</Extensions>
</Filter>
</ItemGroup>
<ItemGroup>
<ClInclude Include="src\opengl\camera.h">
<Filter>Headerdateien</Filter>
</ClInclude>
<ClInclude Include="src\opengl\controller.h">
<Filter>Headerdateien</Filter>
</ClInclude>
<ClInclude Include="src\opengl\shader.h">
<Filter>Headerdateien</Filter>
</ClInclude>
<ClInclude Include="src\opengl\simulation.h">
<Filter>Headerdateien</Filter>
</ClInclude>
<ClInclude Include="src\opengl\visualization.h">
<Filter>Headerdateien</Filter>
</ClInclude>
<ClInclude Include="src\scenarios\SWE_simple_scenarios_vis.h">
<Filter>Headerdateien</Filter>
</ClInclude>
<ClInclude Include="src\scenarios\SWE_VisInfo.h">
<Filter>Headerdateien</Filter>
</ClInclude>
<ClInclude Include="src\scenarios\SWE_VtkScenario.h">
<Filter>Headerdateien</Filter>
</ClInclude>
<ClInclude Include="src\scenarios\SWE_VtkScenarioVisInfo.h">
<Filter>Headerdateien</Filter>
</ClInclude>
</ItemGroup>
<ItemGroup>
<ClCompile Include="src\examples\swe_opengl.cpp">
<Filter>Quelldateien</Filter>
</ClCompile>
<ClCompile Include="src\opengl\camera.cpp">
<Filter>Quelldateien</Filter>
</ClCompile>
<ClCompile Include="src\opengl\controller.cpp">
<Filter>Quelldateien</Filter>
</ClCompile>
<ClCompile Include="src\opengl\shader.cpp">
<Filter>Quelldateien</Filter>
</ClCompile>
<ClCompile Include="src\opengl\visualization.cpp">
<Filter>Quelldateien</Filter>
</ClCompile>
<ClCompile Include="src\scenarios\SWE_VtkScenario.cpp">
<Filter>Quelldateien</Filter>
</ClCompile>
<ClCompile Include="src\scenarios\SWE_VtkScenarioVisInfo.cpp">
<Filter>Quelldateien</Filter>
</ClCompile>
<ClCompile Include="src\SWE_Block.cpp">
<Filter>Quelldateien</Filter>
</ClCompile>
<ClCompile Include="src\SWE_WavePropagationBlock.cpp">
<Filter>Quelldateien</Filter>
</ClCompile>
</ItemGroup>
<ItemGroup>
<CudaCompile Include="src\opengl\simulation.cu">
<Filter>Quelldateien</Filter>
</CudaCompile>
<CudaCompile Include="src\SWE_RusanovBlockCUDA_kernels.cu">
<Filter>Quelldateien</Filter>
</CudaCompile>
<CudaCompile Include="src\SWE_BlockCUDA.cu">
<Filter>Quelldateien</Filter>
</CudaCompile>
<CudaCompile Include="src\SWE_BlockCUDA_kernels.cu">
<Filter>Quelldateien</Filter>
</CudaCompile>
<CudaCompile Include="src\SWE_RusanovBlockCUDA.cu">
<Filter>Quelldateien</Filter>
</CudaCompile>
<CudaCompile Include="src\SWE_WavePropagationBlockCuda.cu">
<Filter>Quelldateien</Filter>
</CudaCompile>
<CudaCompile Include="src\SWE_WavePropagationBlockCuda_kernels.cu">
<Filter>Quelldateien</Filter>
</CudaCompile>
</ItemGroup>
<ItemGroup>
<None Include="src\SWE_RusanovBlockCUDA_kernels.hh">
<Filter>Headerdateien</Filter>
</None>
<None Include="src\SWE_RusanovBlockCUDA.hh">
<Filter>Headerdateien</Filter>
</None>
<None Include="src\SWE_BlockCUDA_kernels.hh">
<Filter>Headerdateien</Filter>
</None>
<None Include="src\SWE_BlockCUDA.hh">
<Filter>Headerdateien</Filter>
</None>
<None Include="src\SWE_Block.hh">
<Filter>Headerdateien</Filter>
</None>
<None Include="src\SWE_WavePropagationBlock.hh">
<Filter>Headerdateien</Filter>
</None>
<None Include="src\SWE_WavePropagationBlockCuda.hh">
<Filter>Headerdateien</Filter>
</None>
<None Include="src\SWE_WavePropagationBlockCuda_kernels.hh">
<Filter>Headerdateien</Filter>
</None>
</ItemGroup>
</Project>
\ No newline at end of file
<?xml version="1.0" encoding="utf-8"?>
<Project ToolsVersion="4.0" xmlns="http://schemas.microsoft.com/developer/msbuild/2003">
</Project>
\ No newline at end of file
......@@ -3,6 +3,10 @@ Microsoft Visual Studio Solution File, Format Version 11.00
# Visual C++ Express 2010
Project("{8BC9CEB8-8B4A-11D0-8D11-00A0C91BC942}") = "SWE for Windows", "SWE for Windows.vcxproj", "{2582913D-07C3-6F6B-46E0-6B3942A726C4}"
EndProject
Project("{8BC9CEB8-8B4A-11D0-8D11-00A0C91BC942}") = "SWE OpenGL", "SWE OpenGL.vcxproj", "{88B845B4-569C-49D1-9C07-967D2E58EAC9}"
EndProject
#Project("{8BC9CEB8-8B4A-11D0-8D11-00A0C91BC942}") = "SWE OpenGL", "..\..\..\visual studio 2010\Projects\SWE #OpenGL\SWE OpenGL.vcxproj", "{88B845B4-569C-49D1-9C07-967D2E58EAC9}"
#EndProject
Global
GlobalSection(SolutionConfigurationPlatforms) = preSolution
CUDA|Win32 = CUDA|Win32
......@@ -16,6 +20,16 @@ Global
{2582913D-07C3-6F6B-46E0-6B3942A726C4}.Debug|Win32.Build.0 = Debug|Win32
{2582913D-07C3-6F6B-46E0-6B3942A726C4}.Release|Win32.ActiveCfg = Release|Win32
{2582913D-07C3-6F6B-46E0-6B3942A726C4}.Release|Win32.Build.0 = Release|Win32
{88B845B4-569C-49D1-9C07-967D2E58EAC9}.Debug|Win32.ActiveCfg = Debug|Win32
{88B845B4-569C-49D1-9C07-967D2E58EAC9}.Debug|Win32.Build.0 = Debug|Win32
{88B845B4-569C-49D1-9C07-967D2E58EAC9}.Release|Win32.ActiveCfg = Release|Win32
{88B845B4-569C-49D1-9C07-967D2E58EAC9}.Release|Win32.Build.0 = Release|Win32
# {88B845B4-569C-49D1-9C07-967D2E58EAC9}.CUDA|Win32.ActiveCfg = Release|Win32
# {88B845B4-569C-49D1-9C07-967D2E58EAC9}.CUDA|Win32.Build.0 = Release|Win32
# {88B845B4-569C-49D1-9C07-967D2E58EAC9}.Debug|Win32.ActiveCfg = Debug|Win32
# {88B845B4-569C-49D1-9C07-967D2E58EAC9}.Debug|Win32.Build.0 = Debug|Win32
# {88B845B4-569C-49D1-9C07-967D2E58EAC9}.Release|Win32.ActiveCfg = Release|Win32
# {88B845B4-569C-49D1-9C07-967D2E58EAC9}.Release|Win32.Build.0 = Release|Win32
EndGlobalSection
GlobalSection(SolutionProperties) = preSolution
HideSolutionNode = FALSE
......
......@@ -102,6 +102,9 @@
<OptimizeReferences>true</OptimizeReferences>
<AdditionalDependencies>cudart.lib;%(AdditionalDependencies)</AdditionalDependencies>
</Link>
<CudaCompile>
<CodeGeneration>compute_20,sm_20</CodeGeneration>
</CudaCompile>
</ItemDefinitionGroup>
<ItemGroup>
<ClCompile Include="src\examples\swe_wavepropagation.cpp" />
......@@ -113,6 +116,12 @@
<ClInclude Include="src\scenarios\SWE_AsagiScenario.hpp" />
<ClInclude Include="src\scenarios\SWE_Scenario.h" />
<ClInclude Include="src\scenarios\SWE_simple_scenarios.h" />
<ClInclude Include="src\solvers\augrie.hpp" />
<ClInclude Include="src\solvers\augriegeoclaw.hpp" />
<ClInclude Include="src\solvers\fwave.hpp" />
<ClInclude Include="src\solvers\FWaveCuda.h" />
<ClInclude Include="src\solvers\hybrid.hpp" />
<ClInclude Include="src\solvers\wavepropagation.hpp" />
<ClInclude Include="src\tools\Logger.hpp" />
</ItemGroup>
<ItemGroup>
......@@ -144,6 +153,15 @@
<FileType>Document</FileType>
</ClInclude>
</ItemGroup>
<ItemGroup>
<None Include="src\swe_block.hh" />
<None Include="src\swe_blockcuda.hh" />
<None Include="src\SWE_BlockCUDA_kernels.hh" />
<None Include="src\swe_rusanovblock.hh" />
<None Include="src\SWE_RusanovBlockCUDA.hh" />
<None Include="src\SWE_RusanovBlockCUDA_kernels.hh" />
<None Include="src\swe_wavepropagationblock.hh" />
</ItemGroup>
<Import Project="$(VCTargetsPath)\Microsoft.Cpp.targets" />
<ImportGroup Label="ExtensionTargets">
<Import Project="$(VCTargetsPath)\BuildCustomizations\CUDA 4.2.targets" />
......
......@@ -43,6 +43,24 @@
</ClInclude>
<ClInclude Include="src\SWE_WavePropagationBlockCuda.hh" />
<ClInclude Include="src\SWE_WavePropagationBlockCuda_kernels.hh" />
<ClInclude Include="src\solvers\augrie.hpp">
<Filter>Header Files</Filter>
</ClInclude>
<ClInclude Include="src\solvers\augriegeoclaw.hpp">
<Filter>Header Files</Filter>
</ClInclude>
<ClInclude Include="src\solvers\fwave.hpp">
<Filter>Header Files</Filter>
</ClInclude>
<ClInclude Include="src\solvers\FWaveCuda.h">
<Filter>Header Files</Filter>
</ClInclude>
<ClInclude Include="src\solvers\hybrid.hpp">
<Filter>Header Files</Filter>
</ClInclude>
<ClInclude Include="src\solvers\wavepropagation.hpp">
<Filter>Header Files</Filter>
</ClInclude>
</ItemGroup>
<ItemGroup>
<CudaCompile Include="src\SWE_BlockCUDA.cu" />
......@@ -52,4 +70,13 @@
<CudaCompile Include="src\SWE_WavePropagationBlockCuda.cu" />
<CudaCompile Include="src\SWE_WavePropagationBlockCuda_kernels.cu" />
</ItemGroup>
<ItemGroup>
<None Include="src\swe_blockcuda.hh" />
<None Include="src\swe_rusanovblock.hh" />
<None Include="src\SWE_BlockCUDA_kernels.hh" />
<None Include="src\swe_block.hh" />
<None Include="src\SWE_RusanovBlockCUDA.hh" />
<None Include="src\SWE_RusanovBlockCUDA_kernels.hh" />
<None Include="src\swe_wavepropagationblock.hh" />
</ItemGroup>
</Project>
\ No newline at end of file
......@@ -4,4 +4,7 @@
<LocalDebuggerCommandArguments>160 160 swe_output</LocalDebuggerCommandArguments>
<DebuggerFlavor>WindowsLocalDebugger</DebuggerFlavor>
</PropertyGroup>
<PropertyGroup>
<ShowAllFiles>true</ShowAllFiles>
</PropertyGroup>
</Project>
\ No newline at end of file
......@@ -122,14 +122,6 @@ SWE_BlockCUDA::SWE_BlockCUDA(float _offsetX, float _offsetY, const int i_cudaDev
cudaMalloc((void**)&bd, size);
checkCUDAError("allocate device memory for bd");
//----> COULD BE IMPLEMENTED AS (INEFFICIENT) REFERENCE IMPLEMENTATION
// // allocate CUDA unknowns: maxmimum height and velocity
// size = (nx/TILE_SIZE)*(ny/TILE_SIZE)*sizeof(float);
// cudaMalloc((void**)&maxhd, size);
// checkCUDAError("allocate device memory for maxhd");
// cudaMalloc((void**)&maxvd, size);
// checkCUDAError("allocate device memory for maxvd");
// allocate consecutive memory for 2 columns with three unknowns each
// (h, hu, hv, excluding b) for copy/ghost layer at bottom/top boundary
size = nx+2;
......
......@@ -205,27 +205,3 @@ void kernelTopCopyLayer(float* hd, float* hud, float* hvd,
}
// //******************************************************************
// // kernels to implement boundary conditions
// //******************************************************************
//
//
// /**
// * CUDA kernel for maximum reduction
// * required to compute maximum water height and velocities to determine
// * allow time step
// */
// __global__
// void kernelMaximum(float* maxhd, float* maxvd, int start, int size) {
// int tx = start+threadIdx.x;
// for (int i=size>>1; i>0; i>>=1) {
// __syncthreads();
// if (tx < i) {
// if( maxhd[tx] < maxhd[tx+i] ) maxhd[tx] = maxhd[tx+i];
// if( maxvd[tx] < maxvd[tx+i] ) maxvd[tx] = maxvd[tx+i];
// };
// };
// }
//
//
......@@ -168,6 +168,35 @@ void SWE_WavePropagationBlockCuda::simulateTimestep(float i_dT) {
updateUnknowns(i_dT);
}
/**
* perform forward-Euler time steps, starting with simulation time tStart,:
* until simulation time tEnd is reached;
* device-global variables hd, hud, hvd are updated;
* unknowns h, hu, hv in main memory are not updated.
* Ghost layers and bathymetry sources are updated between timesteps.
* intended as main simulation loop between two checkpoints
*/
__host__
float SWE_WavePropagationBlockCuda::simulate(float tStart, float tEnd) {
float t = tStart;
do {
// set values in ghost cells:
setGhostLayer();
// Compute the numerical fluxes/net-updates in the wave propagation formulation.
computeNumericalFluxes();
// Update the unknowns with the net-updates.
updateUnknowns(maxTimestep);
t += maxTimestep;
// cout << "Simulation at time " << t << endl << flush;
} while(t < tEnd);
return t;
}
/**
* Compute the numerical fluxes (net-update formulation here) on all edges.
*
......
......@@ -74,11 +74,8 @@ class SWE_WavePropagationBlockCuda: public SWE_BlockCUDA {
// compute a single time step (net-updates + update of the cells).
void simulateTimestep( float i_dT );
// TODO: Not implemented
float simulate(float, float) {
assert(false);
return 0;
};
// simulate multiple time steps (start and end time provided as parameters)
float simulate(float, float);
// TODO: not implemented, max time step reduction is done in each call of computeNumericalFluxes(...)
void computeMaxTimestep() {
......
// =====================================================================
// This file is part of SWE_CUDA (see file SWE_Block.cpp for details).
//
// Copyright (C) 2010,2011 Tobias Schnabel
//
// SWE_CUDA is free software: you can redristribute 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.
//
// SWE_CUDA 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 SWE_CUDA. If not, see <http://www.gnu.org/licenses/>.
// =====================================================================
// Project files
#include "SDL.h"
#include "../opengl/simulation.h"
#include "../opengl/visualization.h"
#include "../opengl/controller.h"
#include "../scenarios/SWE_Scenario.h"
#include "../scenarios/SWE_simple_scenarios_vis.h"
#include "../scenarios/SWE_VtkScenarioVisInfo.h"
#include "../SWE_BlockCUDA.hh"
// #include "../SWE_RusanovBlockCUDA.hh"
#include "../SWE_WavePropagationBlockCUDA.hh"
// For SDL compatibility
#undef main
// Display settings
#define SCREEN_WIDTH 800
#define SCREEN_HEIGHT 600
// Number of nodes (not cells) of grid
#define GRID_XSIZE 401
#define GRID_YSIZE 401
#define WINDOW_TITLE "Shallow Water Equations v1.2"
/**
Main routine.
Takes 1 command line parameter, which specifies optional input file
Control keys:
- [ESC]: quit
- [SPACE]: pause/resume simulation
- [RIGHT]: advance single step (only when paused)
- [s]: save current simulation
- [r]: reset simulation
- [w]: toggle rendering mode
Mouse:
- Left button: rotating
- Right button: panning
- Scroll wheel: zooming
The class concept relies on the Model-View-Controller (MVC)
architecture.
The controller class handles all user input and updates the
visualization and simulation correspondingly.
The simulation (model) advances stepwise, the visualization (view)
displays then the updated data.
*/
int main(int argc, char *argv[])
{
int done = 0;
int gridx = GRID_XSIZE;
int gridy = GRID_YSIZE;
int nx = gridx-1;
int ny = gridy-1;
printf("Starting up...\n");
// Initialize visualization
Visualization visualization(SCREEN_WIDTH, SCREEN_HEIGHT, WINDOW_TITLE, gridx, gridy);
printf("Initialized OpenGL window...\n\n");
// Initialize scenario:
SWE_Scenario* scene = NULL;
SWE_VisInfo* visInfo = NULL;
SWE_BlockCUDA* splash = NULL;
// If input file specified, then read from VTK file:
if (argc > 1) {
SWE_VtkScenarioVisInfo* newScene = SWE_VtkScenarioVisInfo::readVtkFile(argv[1]);
// NOTE: Simulation uses a fixed resolution (independent of VTK file)
scene = newScene;
visInfo = newScene;
printf("Scenario read from input file %s\n\n", argv[1]);
};
if (scene == NULL) {
// ... if VTK file not specified (or was not read successfully)
// use splashing pool scenario ...
SWE_SplashingPoolScenarioVisInfo* newScene = new SWE_SplashingPoolScenarioVisInfo();
scene = newScene;
visInfo = newScene;
};
// define grid size and initial time step
float dx = (scene->getBoundaryPos(BND_RIGHT) - scene->getBoundaryPos(BND_LEFT) )/nx;
float dy = (scene->getBoundaryPos(BND_TOP) - scene->getBoundaryPos(BND_BOTTOM) )/ny;
SWE_Block::initGridData(nx,ny,dx,dy);
// splash = new SWE_RusanovBlockCUDA();
splash = new SWE_WavePropagationBlockCuda();
// define boundary type at all four domain boundaries:
splash->setWallBoundaries(); // walls at all boundaries
// Initialize simulation
printf("Init simulation\n\n");
Simulation sim(nx,ny,dx,dy, scene, visInfo, splash);
printf("Init visualisation\n\n");
visualization.init(&sim);
// Initialize controller
Controller controller(&sim, &visualization);
printf("Start simulation\n\n");
// sim.saveToFile();
// Enter the main loop
while ( ! done )
{
// Handle events
done = controller.handleEvents();
if (controller.isPaused()) {
// We're paused, only render new display
visualization.renderDisplay();
}
else if (controller.hasFocus()) {
// Simulate, update visualization data
sim.runCuda(visualization.getCudaWaterSurfacePtr(), visualization.getCudaNormalsPtr());
// Render new data
visualization.renderDisplay();
}
}
// Clean everything up
visualization.cleanUp();
// delete scene;
// delete splash;
return 0;
}
// =====================================================================
// This file is part of SWE_CUDA (see file SWE_Block.cu for details).
//
// Copyright (C) 2010,2011 Tobias Schnabel
//
// SWE_CUDA is free software: you can redristribute 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.
//
// SWE_CUDA 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 SWE_CUDA. If not, see <http://www.gnu.org/licenses/>.
// =====================================================================
#include "camera.h"
#include <sstream>
/**
Constructor
@param view_distance initial view distance from the origin
@param window_title title of the current window
*/
Camera::Camera(float view_distance, const char* window_title) {
// Initialize member variables
cameraX = -0.3f;
cameraY = 1.0f;
cameraZ = 2.0f;
objectX = 0.0f;
objectY = 0.0f;
objectZ = 0.0f;
angleX = 0.0f;
angleY = 0.0f;
zoomfactor = 0.85f*view_distance;
win_title = window_title;
frames = 0;
lastTime = 0;
oldMouseX = 0;
oldMouseY = 0;
// Reset framebuffer
glClear(GL_COLOR_BUFFER_BIT | GL_DEPTH_BUFFER_BIT);
glMatrixMode(GL_MODELVIEW);
glLoadIdentity();
}
/**
Set the camera via gluLookAt and set the light position afterwards
*/
void Camera::setCamera() {
// Clear our buffers
glClear(GL_COLOR_BUFFER_BIT | GL_DEPTH_BUFFER_BIT);
glMatrixMode(GL_MODELVIEW);
glLoadIdentity();
// Set camera position, center of scene and up vector
gluLookAt(cameraX*zoomfactor, cameraY*zoomfactor, cameraZ*zoomfactor,
0,0,0, // center of scene
0, 1, 0); // up vector
// Translate and rotate our object
glTranslated(objectX, objectY, objectZ);
this->rotateObject();
glScalef(1.0f, 1.0f, -1.0f);
// Position the light source
// Note: Coordinates are multiplied by current modelview-matrix
// by OpenGL
GLfloat LightPosition[]= { 1.0f, 1.0f, 0.0f, 0.0f };
glLightfv(GL_LIGHT0, GL_POSITION,LightPosition);
}
/**
Rotate camera
*/
void Camera::rotateObject() {
glRotatef(angleX, 1.0, 0.0f, 0.0f);
glRotatef(angleY, 0.0f, 1.0f, 0.0f);
}
/**
Increment viewing orientation of the camera
@param angle_dX angle relative to the x-axis
@param angle_dY angle relative to the rotated y-axis
*/
void Camera::orient( float angle_dX, float angle_dY )
{
angleX += angle_dX;
angleY += angle_dY;
}
/**
Zoom in
@param scaleFactor factor which is used for zooming
*/
void Camera::zoomIn( float scaleFactor )
{
if (zoomfactor > 3.0f) {
zoomfactor /= scaleFactor;
}
}
/**
Zoom out
@param scaleFactor factor which is used for zooming
*/
void Camera::zoomOut( float scaleFactor )
{
zoomfactor *= scaleFactor;
}
/**
Calculates the current framerate, updates the window title and
swaps framebuffers to display the new image
*/
void Camera::displayImage() {
// Set up stringstream for window title
std::stringstream ss;
ss.setf(std::ios::fixed, std::ios::floatfield);
ss.setf(std::ios::showpoint);
ss.precision(1);
// Calculate framerate
float fps;
frames++;
if (lastTime == 0) {
lastTime = SDL_GetTicks();
} else if( (SDL_GetTicks() - lastTime) >= 500) {
fps = (float)frames/(SDL_GetTicks() - lastTime)*1000.0f;
ss << win_title << " (" << fps << " fps)";
SDL_WM_SetCaption(ss.str().c_str(), NULL);
frames = 0;
lastTime = SDL_GetTicks();
}
// Swap frame buffer
SDL_GL_SwapBuffers();
}
/**
User starts dragging.
Remember the old mouse coordinates.
*/
void Camera::startPanning(int xPos, int yPos) {
oldMouseX = xPos;
oldMouseY = yPos;
}
/**
User drags our object.
Transform screen coordinates into world coordinates
and update the objects position
*/
void Camera::panning(int newX, int newY) {
GLint viewport[4];
GLdouble modelview[16];
GLdouble projection[16];
GLfloat winX, winY, winZ;
GLdouble posX, posY, posZ;
GLdouble pos2X, pos2Y, pos2Z;
GLdouble dX, dY, dZ;
// Draw invisible fake plane
setCamera();
glColorMask( GL_FALSE, GL_FALSE, GL_FALSE, GL_FALSE );
GLfloat dim = 100000.0f;
glBegin(GL_QUADS);
glColor3f(1.0f, 1.0f, 0.3f);
glVertex3f( dim,0.0f, dim);
glVertex3f(-dim,0.0f, dim);
glVertex3f(-dim,0.0f, -dim);
glVertex3f( dim,0.0f, -dim);
glEnd();
glColorMask(GL_TRUE,GL_TRUE,GL_TRUE,GL_TRUE);
// Get projection matrices
glGetDoublev( GL_MODELVIEW_MATRIX, modelview );
glGetDoublev( GL_PROJECTION_MATRIX, projection );
glGetIntegerv( GL_VIEWPORT, viewport );
// Get drag positions in object coordinates
winX = (GLfloat)newX;
winY = (GLfloat)viewport[3] - (GLfloat)newY;
glReadPixels( newX, GLint(winY), 1, 1, GL_DEPTH_COMPONENT, GL_FLOAT, &winZ );
gluUnProject( winX, winY, winZ, modelview, projection, viewport, &posX, &posY, &posZ);
winX = (GLfloat)oldMouseX;
winY = (GLfloat)viewport[3] - (GLfloat)oldMouseY;
glReadPixels( oldMouseX, GLint(winY), 1, 1, GL_DEPTH_COMPONENT, GL_FLOAT, &winZ );
gluUnProject( winX, winY, winZ, modelview, projection, viewport, &pos2X, &pos2Y, &pos2Z);
glClear(GL_DEPTH_BUFFER_BIT);
// Calculate drag distances
dX = (posX-pos2X);
dY = (posY-pos2Y);
dZ = (posZ-pos2Z);
// Convert drag distances to world coordinates
glPushMatrix();
glLoadIdentity();
rotateObject();
glScalef(1.0f, 1.0f, -1.0f);
glGetDoublev( GL_MODELVIEW_MATRIX, modelview );
GLdouble x = modelview[0]*dX + modelview[4]*dY + modelview[8]*dZ;
GLdouble y = modelview[1]*dX + modelview[5]*dY + modelview[9]*dZ;
GLdouble z = modelview[2]*dX + modelview[6]*dY + modelview[10]*dZ;
glPopMatrix();
// Set the new world coordinates of our object
objectX += x;
objectY += y;
objectZ += z;
// Save our current mouse position
oldMouseX = newX;
oldMouseY = newY;
}
\ No newline at end of file
#ifndef CAMERA_H
#define CAMERA_H
// =====================================================================
// This file is part of SWE_CUDA (see file SWE_Block.cu for details).
//
// Copyright (C) 2010,2011 Tobias Schnabel
//
// SWE_CUDA is free software: you can redristribute 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.
//
// SWE_CUDA 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 SWE_CUDA. If not, see <http://www.gnu.org/licenses/>.
// =====================================================================
#include <SDL.h>
#include <SDL_opengl.h>
class Camera {
public:
Camera(float view_distance, const char* window_title);
// Set modelview matrix
void setCamera ();
// Change viewing properties
void orient( float angX, float angY );
void zoomIn( float scaleFactor );
void zoomOut( float scaleFactor );
void startPanning(int xPos, int yPos);
void panning(int newX, int newY);
// Update framebuffer
void displayImage();
private:
// Position of the camera
float cameraX;
float cameraY;
float cameraZ;
// Position of the object
GLdouble objectX;
GLdouble objectY;
GLdouble objectZ;
// Zoom factor
float zoomfactor;
float angleX, angleY;
// Helper variables
unsigned int frames;
unsigned int lastTime;
unsigned int oldMouseX, oldMouseY;
unsigned int newMouseX, newMouseY;
// Window title
const char* win_title;
void rotateObject();
};
#endif
\ No newline at end of file
// =====================================================================
// This file is part of SWE_CUDA (see file SWE_Block.cu for details).
//
// Copyright (C) 2010,2011 Tobias Schnabel
//
// SWE_CUDA is free software: you can redristribute 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.
//
// SWE_CUDA 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 SWE_CUDA. If not, see <http://www.gnu.org/licenses/>.
// =====================================================================
#include "controller.h"
#include "../scenarios/SWE_simple_scenarios_vis.h"
/**
Constructor
@param sim instance of simulation class
@param vis instance of visualization class
*/
Controller::Controller(Simulation* sim, Visualization* vis) {
simulation = sim;
visualization = vis;
isActive = true;
done = false;
paused = false;
allowStep = false;
}
/**
Process all user events in a loop
Returns true, when user wants to quit
*/
bool Controller::handleEvents() {
// Handle events
SDL_Event event;
allowStep = false;
while ( SDL_PollEvent(&event) )
{
switch( event.type )
{
case SDL_MOUSEMOTION:
if (event.motion.state & SDL_BUTTON(SDL_BUTTON_LEFT)) {
// Left mouse button dragging
visualization->camera->orient(0.5f*event.motion.yrel, 0.5f*event.motion.xrel);
}
else if (event.motion.state & SDL_BUTTON(SDL_BUTTON_RIGHT)) {
// Panning happens here
visualization->camera->panning(event.motion.x, event.motion.y);
}
break;
case SDL_MOUSEBUTTONDOWN:
// Scroll wheel down
if ( event.button.button == SDL_BUTTON_WHEELDOWN ) {
// Zoom out
visualization->camera->zoomOut(1.2f);
} else if (event.button.button == SDL_BUTTON_RIGHT) {
visualization->camera->startPanning(event.button.x, event.button.y);
}
break;
case SDL_MOUSEBUTTONUP:
// Scroll wheel up
if ( event.button.button == SDL_BUTTON_WHEELUP ) {
// Zoom in
visualization->camera->zoomIn(1.15f);
}
break;
case SDL_ACTIVEEVENT:
// Don't draw anything if we're getting minimized
if ( event.active.state & SDL_APPACTIVE ) {
isActive = (event.active.gain != 0);
}
break;
case SDL_VIDEORESIZE:
// Our window was resized
visualization->resizeWindow(event.resize.w, event.resize.h);
break;
case SDL_KEYDOWN:
// User has pressed a key
done = handleKeyPress( &event.key.keysym);
break;
case SDL_QUIT:
// Window close request
done = true;
break;
default:
break;
}
}
return done;
}
/**
Returns true, when window has focus
*/
bool Controller::hasFocus() {
return isActive && !done;
}
/**
Return whether program is currently paused
*/
bool Controller::isPaused() {
return paused && !allowStep;
}
/**
Process single keyboard event
@param *keysym pointer to the sdl keyevent structure
*/
bool Controller::handleKeyPress( SDL_keysym *keysym) {
switch ( keysym->sym ) {
case SDLK_ESCAPE:
// ESC key was pressed -> exit
return true;
break;
case SDLK_r:
// Restart simulation
allowStep = paused;
simulation->restart();
break;
case SDLK_w:
// Switch rendering mode
visualization->toggleRenderingMode();
break;
case SDLK_s:
// Save simulation data to file
simulation->saveToFile();
break;
case SDLK_RIGHT:
// Advance single timestep when paused
allowStep = paused;
break;
case SDLK_SPACE:
// Pause/Resume
paused = !paused;
break;
case SDLK_1:
// Load scenario 1
{
allowStep = paused;
SWE_RadialDamBreakScenarioVisInfo*
newScene = new SWE_RadialDamBreakScenarioVisInfo();
SWE_VisInfo* visInfo = newScene;
// define grid size and initial time step
float dx = (newScene->getBoundaryPos(BND_RIGHT) - newScene->getBoundaryPos(BND_LEFT) )/SWE_Block::getNx();
float dy = (newScene->getBoundaryPos(BND_TOP) - newScene->getBoundaryPos(BND_BOTTOM) )/SWE_Block::getNy();
SWE_Block::initGridData(SWE_Block::getNx(),SWE_Block::getNy(),dx,dy);
simulation->loadNewScenario(newScene, visInfo);
visualization->updateBathymetryVBO(simulation);
}
break;
case SDLK_2:
// Load scenario 2
{
allowStep = paused;
SWE_Scenario* newScene = new SWE_BathymetryDamBreakScenario;
// define grid size and initial time step
float dx = (newScene->getBoundaryPos(BND_RIGHT) - newScene->getBoundaryPos(BND_LEFT) )/SWE_Block::getNx();
float dy = (newScene->getBoundaryPos(BND_TOP) - newScene->getBoundaryPos(BND_BOTTOM) )/SWE_Block::getNy();
SWE_Block::initGridData(SWE_Block::getNx(),SWE_Block::getNy(),dx,dy);
simulation->loadNewScenario(newScene, NULL);
visualization->updateBathymetryVBO(simulation);
}
break;
case SDLK_3:
// Load scenario 3
{
allowStep = paused;
SWE_SplashingPoolScenarioVisInfo*
newScene = new SWE_SplashingPoolScenarioVisInfo;
// define grid size and initial time step
float dx = (newScene->getBoundaryPos(BND_RIGHT) - newScene->getBoundaryPos(BND_LEFT) )/SWE_Block::getNx();
float dy = (newScene->getBoundaryPos(BND_TOP) - newScene->getBoundaryPos(BND_BOTTOM) )/SWE_Block::getNy();
SWE_Block::initGridData(SWE_Block::getNx(),SWE_Block::getNy(),dx,dy);
simulation->loadNewScenario(newScene, newScene);
visualization->updateBathymetryVBO(simulation);
}
break;
default:
break;
}
return false;
}
#ifndef CONTROLLER_H
#define CONTROLLER_H
// =====================================================================
// This file is part of SWE_CUDA (see file SWE_Block.cu for details).
//
// Copyright (C) 2010,2011 Tobias Schnabel
//
// SWE_CUDA is free software: you can redristribute 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.
//
// SWE_CUDA 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 SWE_CUDA. If not, see <http://www.gnu.org/licenses/>.
// =====================================================================
#include <SDL.h>
#include "simulation.h"
#include "visualization.h"
class Controller {
public:
Controller(Simulation* sim, Visualization* vis);
// Process new events
bool handleEvents();
// Return if window is enabled
bool hasFocus();
// Return if program is paused
bool isPaused();
private:
// Internal functions
bool isActive;
bool done;
bool paused;
bool allowStep;
// References to other classes
Simulation* simulation;
Visualization* visualization;
// Handle keyboard events
bool handleKeyPress( SDL_keysym *keysym);
};
#endif
\ No newline at end of file
// =====================================================================
// This file is part of SWE_CUDA (see file SWE_Block.cu for details).
//
// Copyright (C) 2010,2011 Tobias Schnabel
//
// SWE_CUDA is free software: you can redristribute 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.
//
// SWE_CUDA 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 SWE_CUDA. If not, see <http://www.gnu.org/licenses/>.
// =====================================================================
varying vec3 V;
varying vec3 N;
varying vec4 ambient;
/* Shading of water surfaces with fresnel effect */
void main()
{
// Helper variable for color
vec4 color = ambient;
// Get normal and view vector
vec3 view = normalize(V);
vec3 normal = normalize(N); ;
if (!gl_FrontFacing) {
// Handle backface rendering
normal = -normal;
}
// Compute fresnel factor
float fresnelFactor = 0.1+0.1*pow(1+dot(normal,view), 3.0);
// Combine color with fresnel factor
color += gl_LightSource[0].specular * fresnelFactor;
color = mix(color, vec4(0.2, 0.4, 0.9, 1.0), 0.2);
// Set constant transparency
gl_FragColor = vec4(color.xyz, 0.7);
}
// =====================================================================
// This file is part of SWE_CUDA (see file SWE_Block.cu for details).
//
// Copyright (C) 2010,2011 Tobias Schnabel
//
// SWE_CUDA is free software: you can redristribute 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.
//
// SWE_CUDA 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 SWE_CUDA. If not, see <http://www.gnu.org/licenses/>.
// =====================================================================
#include "shader.h"
/**
Constructor.
Check whether shaders are supported. If yes, load vertex and fragment
shader from textfile into memory and compile
@param vertexShaderFile name of the text file containing the vertex
shader code
@param fragmentShaderFile name of the text file containing the fragment
shader code
*/
Shader::Shader(char const * vertexShaderFile, char const * fragmentShaderFile) {
shdrSupport = isExtensionSupported("GL_ARB_vertex_shader")
&& isExtensionSupported("GL_ARB_shader_objects")
&& isExtensionSupported("GL_ARB_fragment_shader");
shdrLoaded = false;
if (shdrSupport) {
// Load extensions
glCreateShader = (PFNGLCREATESHADERPROC) SDL_GL_GetProcAddress("glCreateShader");
glCreateProgram = (PFNGLCREATEPROGRAMPROC) SDL_GL_GetProcAddress("glCreateProgram");
glAttachShader = (PFNGLATTACHSHADERPROC)SDL_GL_GetProcAddress("glAttachShader");
glCompileShader = (PFNGLCOMPILESHADERPROC)SDL_GL_GetProcAddress("glCompileShader");
glUseProgram = (PFNGLUSEPROGRAMPROC)SDL_GL_GetProcAddress("glUseProgram");
glDetachShader = (PFNGLDETACHSHADERPROC)SDL_GL_GetProcAddress("glDetachShader");
glDeleteShader = (PFNGLDELETESHADERPROC) SDL_GL_GetProcAddress("glDeleteShader");
glLinkProgram = (PFNGLLINKPROGRAMPROC)SDL_GL_GetProcAddress("glLinkProgram");
glShaderSource = (PFNGLSHADERSOURCEPROC)SDL_GL_GetProcAddress("glShaderSource");
glDeleteProgram = (PFNGLDELETEPROGRAMPROC)SDL_GL_GetProcAddress("glDeleteProgram");
glGetObjectParameterivARB = (PFNGLGETOBJECTPARAMETERIVARBPROC) SDL_GL_GetProcAddress("glGetObjectParameterivARB");
glGetShaderiv = (PFNGLGETSHADERIVPROC) SDL_GL_GetProcAddress("glGetShaderiv");
glGetShaderInfoLog = (PFNGLGETSHADERINFOLOGPROC) SDL_GL_GetProcAddress("glGetShaderInfoLog");
glGetProgramInfoLog = (PFNGLGETPROGRAMINFOLOGPROC) SDL_GL_GetProcAddress("glGetProgramInfoLog");
// Read shader files
bool readSuccess = false;
vertexShaderSource = NULL;
fragmentShaderSource = NULL;
readSuccess = readShaderFile(vertexShaderFile, vertexShaderSource, vertexShaderLength);
readSuccess &= readShaderFile(fragmentShaderFile, fragmentShaderSource, fragmentShaderLength);
// Check if reading was completed
if (readSuccess) {
// Create Vertex-Shader
vertexShader = glCreateShader(GL_VERTEX_SHADER);
glShaderSource(vertexShader, 1, const_cast<GLchar const **>(&vertexShaderSource), &vertexShaderLength);
glCompileShader(vertexShader);
// Create Fragment-Shader
fragmentShader = glCreateShader(GL_FRAGMENT_SHADER);
glShaderSource(fragmentShader, 1, const_cast<GLchar const **>(&fragmentShaderSource), &fragmentShaderLength);
glCompileShader(fragmentShader);
if (isShaderCompiled(vertexShader, "Vertex-Shader")
&& isShaderCompiled(fragmentShader, "Fragment-Shader"))
{
// Create Shader program
program = glCreateProgram();
glAttachShader(program, vertexShader);
glAttachShader(program, fragmentShader);
glLinkProgram(program);
shdrLoaded = true;
} else {
// Errors while compiling shaders
shdrSupport = false;
glDeleteShader(vertexShader);
delete[] vertexShaderSource;
glDeleteShader(fragmentShader);
delete[] fragmentShaderSource;
}
} else {
if (vertexShaderSource != NULL) delete[] vertexShaderSource;
if (fragmentShaderSource != NULL) delete[] fragmentShaderSource;
}
}
}
/**
Destructor.
Unload shaders and free resources.
*/
Shader::~Shader() {
// Unload Shaders
if (shdrLoaded) {
glUseProgram(0);
glDeleteProgram(program);
glDetachShader(program, vertexShader);
glDeleteShader(vertexShader);
delete[] vertexShaderSource;
glDetachShader(program, fragmentShader);
glDeleteShader(fragmentShader);
delete[] fragmentShaderSource;
}
}
/**
Replaces OpenGL shaders by our custom shaders
*/
void Shader::enableShader() {
if (shdrLoaded) {
glUseProgram(program);
}
}
/**
Restores OpenGL default shaders
*/
void Shader::disableShader() {
if (shdrLoaded) {
glUseProgram(0);
}
}
/**
Returns, whether shaders are supported by graphics hardware
*/
bool Shader::shadersSupported() {
return shdrSupport;
}
/**
Returns, whether shaders could by loaded successfully
*/
bool Shader::shadersLoaded() {
return shdrLoaded;
}
/**
Returns, whether a special extension is supported by the current
graphics card
@param szTargetExtention string describing the extension to look for
*/
bool Shader::isExtensionSupported(const char* szTargetExtension )
{
const unsigned char *pszExtensions = NULL;
const unsigned char *pszStart;
unsigned char *pszWhere, *pszTerminator;
// Extension names should not have spaces
pszWhere = (unsigned char *) strchr( szTargetExtension, ' ' );
if( pszWhere || *szTargetExtension == '\0' )
return false;
// Get Extensions String
pszExtensions = glGetString( GL_EXTENSIONS );
// Search The Extensions String For An Exact Copy
pszStart = pszExtensions;
for(;;)
{
pszWhere = (unsigned char *) strstr( (const char *) pszStart, szTargetExtension );
if( !pszWhere )
break;
pszTerminator = pszWhere + strlen( szTargetExtension );
if( pszWhere == pszStart || *( pszWhere - 1 ) == ' ' )
if( *pszTerminator == ' ' || *pszTerminator == '\0' )
return true;
pszStart = pszTerminator;
}
return false;
}
/**
Returns, whether a special extension is supported by the current
graphics card
@param filename shader file
@return shaderSource file content
@return length length of the file
@returns bool true, if file has been read successfully
*/
bool Shader::readShaderFile(char const * filename, GLchar * & shaderSource, GLint & length)
{
using namespace std;
string fullFileName = "";
// On unix: determine directory of exec
#ifdef __unix__
// Code taken from:
// http://www.gamedev.net/community/forums/topic.asp?topic_id=459511
std::string path = "";
pid_t pid = getpid();
char buf[20] = {0};
sprintf(buf,"%d",pid);
std::string _link = "/proc/";
_link.append( buf );
_link.append( "/exe");
char proc[512];
int ch = readlink(_link.c_str(),proc,512);
if (ch != -1) {
proc[ch] = 0;
path = proc;
std::string::size_type t = path.find_last_of("/");
path = path.substr(0,t);
}
fullFileName = path + string("/");
#endif
fullFileName.append(filename);
ifstream file(fullFileName.c_str(), ios::in);
length = 0;
if (!file.good())
{
return false;
}
file.seekg(0, ios::end);
length = file.tellg();
file.seekg(ios::beg);
if (length == 0 )
{
return false;
}
shaderSource = new GLchar[length];
if (shaderSource == 0)
{
return false;
}
memset(shaderSource, 0, length);
file.read(shaderSource, length);
file.close();
return true;
}
/**
Returns, whether a shader has been compiled with success
@param shader shader id
@param prefix custom error message
*/
bool Shader::isShaderCompiled(GLuint shader, char const * prefix)
{
using namespace std;
GLint compiled(0);
glGetObjectParameterivARB(shader, GL_COMPILE_STATUS, &compiled);
if (!compiled)
{
GLint maxLength(0);
glGetShaderiv(shader, GL_INFO_LOG_LENGTH , &maxLength);
cerr << prefix
<< " -- Compiler Log: "
<< endl
;
if (maxLength > 1)
{
GLchar * log = new GLchar[maxLength];
glGetShaderInfoLog(shader, maxLength, 0, log);
cerr << log
;
delete[] log;
}
}
return (compiled != 0);
}
/**
Returns, whether a shader has been linked with success
@param shader shader id
@param prefix custom error message
*/
bool Shader::isProgramLinked(GLuint program, char const * prefix)
{
using namespace std;
GLint linked(0);
glGetObjectParameterivARB(program, GL_LINK_STATUS, &linked);
if (!linked)
{
GLint maxLength(0);
glGetShaderiv(program, GL_INFO_LOG_LENGTH , &maxLength);
cerr << prefix
<< " -- Linker Log: "
<< endl
;
if (maxLength > 1)
{
GLchar * log = new GLchar[maxLength];
glGetProgramInfoLog(program, maxLength, 0, log);
cerr << log
;
delete[] log;
}
}
return (linked != 0);
}
#ifndef SHADER_H
#define SHADER_H
// =====================================================================
// This file is part of SWE_CUDA (see file SWE_Block.cu for details).
//
// Copyright (C) 2010,2011 Tobias Schnabel
//
// SWE_CUDA is free software: you can redristribute 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.
//
// SWE_CUDA 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 SWE_CUDA. If not, see <http://www.gnu.org/licenses/>.
// =====================================================================
#include <SDL.h>
#include <SDL_opengl.h>
#include <iostream>
#include <fstream>
class Shader {
public:
// Constructor and Destructor
// Load vertex shader file and fragment shader file
Shader(char const * vertexShaderFile, char const * fragmentShaderFile) ;
~Shader();
// Check if shaders are supported
bool shadersSupported();
bool shadersLoaded();
// Shader control
void enableShader();
void disableShader();
private:
// State flags
bool shdrSupport;
bool shdrLoaded;
// Helper functions
bool readShaderFile(char const * filename, GLchar * & shaderSource, GLint & length);
bool isExtensionSupported(const char* szTargetExtension );
bool isShaderCompiled(GLuint shader, char const * prefix);
bool isProgramLinked(GLuint program, char const * prefix);
// Vertex-Shader
GLuint vertexShader;
GLchar * vertexShaderSource;
GLint vertexShaderLength;
// Fragment-Shader
GLuint fragmentShader;
GLchar * fragmentShaderSource;
GLint fragmentShaderLength;
// Shaders id
GLuint program;
// Shader extension function pointers
PFNGLCREATESHADERPROC glCreateShader;
PFNGLCREATEPROGRAMPROC glCreateProgram;
PFNGLATTACHSHADERPROC glAttachShader;
PFNGLCOMPILESHADERPROC glCompileShader;
PFNGLUSEPROGRAMPROC glUseProgram;
PFNGLDETACHSHADERPROC glDetachShader;
PFNGLDELETESHADERPROC glDeleteShader;
PFNGLLINKPROGRAMPROC glLinkProgram;
PFNGLSHADERSOURCEPROC glShaderSource;
PFNGLDELETEPROGRAMPROC glDeleteProgram;
// Shader objects extension pointers
PFNGLGETOBJECTPARAMETERIVARBPROC glGetObjectParameterivARB;
PFNGLGETSHADERIVPROC glGetShaderiv;
PFNGLGETSHADERINFOLOGPROC glGetShaderInfoLog;
PFNGLGETPROGRAMINFOLOGPROC glGetProgramInfoLog;
};
#endif
// =====================================================================
// This file is part of SWE_CUDA (see file SWE_Block.cu for details).
//
// Copyright (C) 2010,2011 Michael Bader, Kaveh Rahnema, Tobias Schnabel
//
// SWE_CUDA is free software: you can redristribute 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.
//
// SWE_CUDA 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 SWE_CUDA. If not, see <http://www.gnu.org/licenses/>.
// =====================================================================
#include "simulation.h"
#include <stdlib.h>
#include "../SWE_BlockCUDA.hh"
/**
Constructor.
Initializes SWE_BlockCUDA and creates a new instance of it.
@param nx number of grid cells in x-direction
@param ny number of grid cells in y-direction
@param dx size of one cell in x-direction
@param dy size of one cell in y-direction
@param scene instance of SWE_Scenario describing the boundaries
*/
Simulation::Simulation (int nx, int ny, float dx, float dy,
SWE_Scenario* scene, SWE_VisInfo* visInfo,
SWE_BlockCUDA* _splash)
: myScenario(scene),
maxDim( (nx > ny) ? nx : ny ),
maxCellSize( (dx > dy) ? dx : dy ),
splash(_splash),
curTime(0.0f),
isFirstStep(1),
fileNumber(0),
useFileInput(false)
{
Simulation::initBoundaries(myScenario);
if (visInfo != NULL) {
wAverage = visInfo->waterHeightAtRest();
wScale = visInfo->waterVerticalScaling()*(maxDim/20.0f);
wOffset = visInfo->waterDistanceFromGround()*(maxDim/20.0f);
bAverage = visInfo->bathyVerticalCenter();
bScale = visInfo->bathyVerticalScaling()*(maxDim/20.0f);
bOffset = visInfo->bathyDistanceFromGround()*(maxDim/20.0f);
} else {
getScalingApproximation(splash->getWaterHeight(), splash->getBathymetry());
}
}
/**
Destructor.
Delete current SWE_BlockCUDA and SWE_Scenario instance.
*/
Simulation::~Simulation () {
}
void Simulation::loadNewScenario(SWE_Scenario* scene, SWE_VisInfo* visInfo) {
delete myScenario;
myScenario = scene;
curTime = 0.0f;
isFirstStep = 1;
useFileInput = false;
Simulation::initBoundaries(myScenario);
if (visInfo != NULL) {
wAverage = visInfo->waterHeightAtRest();
wScale = visInfo->waterVerticalScaling()*(maxDim/20.0f);
wOffset = visInfo->waterDistanceFromGround()*(maxDim/20.0f);
bAverage = visInfo->bathyVerticalCenter();
bScale = visInfo->bathyVerticalScaling()*(maxDim/20.0f);
bOffset = visInfo->bathyDistanceFromGround()*(maxDim/20.0f);
} else {
getScalingApproximation(splash->getWaterHeight(), splash->getBathymetry());
}
}
/**
This is the main simulation procedure. Simulates one timestep and computes
normals afterwards.
@param vbo_resource cuda resource holding the vertex positions
@param vbo_normals cuda resource holding the normals
*/
void Simulation::runCuda(struct cudaGraphicsResource **vbo_resource, struct cudaGraphicsResource **vbo_normals)
{
// map OpenGL buffer object for writing from CUDA
float3 *dptr, *dptr2;
size_t num_bytes, num_bytes2;
cudaGraphicsMapResources(1, vbo_resource, 0);
cudaGraphicsResourceGetMappedPointer((void **)&dptr, &num_bytes,
*vbo_resource);
// run simulation and fill our vbo
calculateWaterSurface(dptr);
// unmap buffer object
cudaGraphicsUnmapResources(1, vbo_resource, 0);
checkCUDAError("Fehler bei Simulations-VBO\n");
//Calculate normals
cudaGraphicsMapResources(1, vbo_resource, 0);
cudaGraphicsResourceGetMappedPointer((void **)&dptr, &num_bytes,
*vbo_resource);
cudaGraphicsMapResources(1, vbo_normals, 0);
cudaGraphicsResourceGetMappedPointer((void **)&dptr2, &num_bytes2,
*vbo_normals);
// calculate normals of our new water surface
calculateNormals(dptr, dptr2);
// unmap buffer objects
cudaGraphicsUnmapResources(1, vbo_resource, 0);
checkCUDAError("Fehler bei Simulations-VBO\n");
cudaGraphicsUnmapResources(1, vbo_normals, 0);
checkCUDAError("Fehler bei Normalen-VBO\n");
}
/**
Init new boundaries defined by the scene
@param scene instance of SWE_Scenario
*/
void Simulation::initBoundaries(SWE_Scenario* scene) {
std::cout << "Init Scenario\n" << flush;
splash->initScenario(*scene);
splash->setGhostLayer();
}
/**
Restarts the simulation. Restores the initial bondaries.
*/
void Simulation::restart() {
curTime = 0.0f;
isFirstStep = 1;
initBoundaries(myScenario);
}
/**
Sets the bathymetry buffer. Buffer contains vertex position and
vertex normal in sequence.
@param bath float array in which computed values will be stored
*/
void Simulation::setBathBuffer(float* bath) {
// splash->setBathBuffer(output, bAverage, bScale, bOffset);
const Float2D& b = splash->getBathymetry();
int nx = b.getRows()-2;
int ny = b.getCols()-2;
// Set vertex coordinates
for (int j=0; j<ny+1;j++) {
for (int i=0;i<nx+1;i++) {
bath[(j*(ny+1) + i)*6] = (float) i;
bath[(j*(ny+1) + i)*6 + 1]= scaleFunction(0.25f*(b[i][j]+b[i+1][j]+b[i][j+1]+b[i+1][j+1]),
bAverage, bScale, bOffset);
bath[(j*(ny+1) + i)*6 + 2] = (float) j;
bath[(j*(ny+1) + i)*6 + 3] = 0.0f;
bath[(j*(ny+1) + i)*6 + 4] = 0.0f;
bath[(j*(ny+1) + i)*6 + 5] = 0.0f;
}
}
// Calculate normals
for(int j=0; j < ny; j++)
{
for(int i=0; i < nx; i++)
{
// Calculate normal vectors for each triangle
float normal1[3];
float normal2[3];
calculateNormal(&bath[(j*(ny+1) + i)*6],
&bath[((j+1)*(ny+1) + i + 1)*6],
&bath[((j+1)*(ny+1) + i)*6],
normal1);
calculateNormal(&bath[(j*(ny+1) + i)*6],
&bath[(j*(ny+1) + i + 1)*6],
&bath[((j+1)*(ny+1) + i + 1)*6],
normal2);
// Copy normals to array
for (int k=0; k < 3; k++) {
bath[(j*(ny+1) + i)*6 + 3 + k] = (normal1[k]+normal2[k])*0.5f;
}
}
}
// Fill boundary regions
for(int x=0; x < ny; x++) {
for (int i=0; i < 3; i++) {
bath[(ny*(ny+1) + x)*6 + 3 + i] = bath[((ny-1)*(ny+1) + x)*6 + 3 + i];
}
}
for(int y=0; y < nx; y++) {
for (int i=0; i < 3; i++) {
bath[(y*(ny+1) + nx)*6 + 3 + i] = bath[(y*(ny+1) + nx - 1)*6 + 3 + i];
}
}
for (int i=0; i < 3; i++) {
bath[(ny*(ny+1) + nx)*6 + 3 + i] = bath[((ny-1)*(ny+1) + nx - 1)*6 + 3 + i];
}
}
/**
Used for debugging purposes only. Write the current simulation and
visualization buffer to a file.
*/
void Simulation::writeDebugOutput(float3* destBuffer) {
splash->writeVTKFile3D(generateFileName("debugoutput",0));
if (destBuffer != NULL) {
debugVisBuffer(destBuffer);
}
}
/**
Save the current simulation state to a vtk file
*/
void Simulation::saveToFile() {
splash->writeVTKFile3D(generateFileName("simdata", fileNumber));
std::cout << "Saved data to: " << generateFileName("simdata", fileNumber) << std::endl;
fileNumber++;
}
/**
Advance single step in simulation on graphics card
and copy results to visualization buffer
@param destBuffer vertex buffer object holding all water surface data
*/
void Simulation::calculateWaterSurface(float3* destBuffer) {
if (isFirstStep == 1) {
isFirstStep = 0;
} else {
//splash->simulateConstTimestep();
// curTime = splash->simulate(curTime, curTime + 5*splash->getMaxTimestep());
float dt = splash->getMaxTimestep();
splash->simulateTimestep(dt);
curTime += dt;
}
// splash->updateVisBuffer(destBuffer, wAverage, wScale, wOffset);
updateVisBuffer(destBuffer);
}
/**
Computes a first approximation of the scaling values needed
for visualization.
Gets called before simulation starts and determines the average,
mininimum and maximum values of the bathymetry and water surface data.
Uses latter values to estimate the scaling factors.
@param h Fload2D array that holds water height
@param b Fload2D array that holds bathymetry data
*/
void Simulation::getScalingApproximation(const Float2D& h, const Float2D& b) {
// Averages of B and (B + H)
float avgB, avgH;
// Minimum values
float minB, minH;
// Maximum values
float maxB, maxH;
int nx = h.getRows()-2;
int ny = h.getCols()-2;
int maxDim = (nx > ny) ? nx : ny;
avgH = 0.0f;
avgB = 0.0f;
minB = b[1][1];
minH = h[1][1];
maxB = b[1][1];
maxH = h[1][1];
for(int i=1; i<=nx; i++) {
for(int j=1; j<=ny; j++) {
// Update averages
avgH += (h[i][j] + b[i][j])/(nx*ny);
avgB += b[i][j]/(nx*ny);
// Update minima
if ((h[i][j] + b[i][j]) < minH)
minH = (h[i][j] + b[i][j]);
if (b[i][j] < minB)
minB = b[i][j];
// Update maxima
if ((h[i][j] + b[i][j]) > maxH)
maxH = (h[i][j] + b[i][j]);
if (b[i][j] > maxB)
maxB = b[i][j];
}
}
// cout << "max and min bathymetry: " << maxB << " ; " << minB << endl;
// bAverage = avgB;
bAverage = (maxB+minB)/2;
cout << "Average bathymetry: " << bAverage << endl;
if (abs(maxB - minB) > 0.0001f) {
bScale = (maxDim/20.0f)/(maxB - minB);
} else {
bScale = (maxDim/20.0f);
}
cout << "Scaling of bathymetry: " << bScale << endl;
bOffset = bScale*(avgB - minB) + maxDim/15.0f;
// bOffset = minB+100;
cout << "bathymetry offset: " << bOffset << endl;
wAverage = avgH;
cout << "Average water height: " << avgH << endl;
if ((maxH - minH) < 0.0001f) {
wScale = 1.0f/(maxH- minH);
} else {
wScale = 1.0f;
}
wScale = (maxDim/40.0)*wScale;
cout << "Scaling of water level: " << wScale << endl;
wOffset = bOffset + (maxDim/5.0f);
cout << "water level offset: " << wOffset << endl;
}
/**
Compute normal of a triangle
@param fVert1-fVert3 vertices of the triangle
@param fNormal resulting normal
*/
void Simulation::calculateNormal(float fVert1[], float fVert2[],
float fVert3[], float fNormal[]) {
float Qx, Qy, Qz, Px, Py, Pz;
Qx = fVert2[0]-fVert1[0];
Qy = fVert2[1]-fVert1[1];
Qz = fVert2[2]-fVert1[2];
Px = fVert3[0]-fVert1[0];
Py = fVert3[1]-fVert1[1];
Pz = fVert3[2]-fVert1[2];
fNormal[0] = Py*Qz - Pz*Qy;
fNormal[1] = Pz*Qx - Px*Qz;
fNormal[2] = Px*Qy - Py*Qx;
}
//==================================================================
// member functions for visualisation and output
//==================================================================
/**
Scale function used for visualization
*/
__device__ __host__
float scaleFunction(float val, float average, float scale, float offset) {
return (val - average)*scale + offset;
}
// /**
// Return coordinates of visBuffer[x][y]
// @param x,y coordinates
// @param ny length in y-direction
// */
// __device__
// int getVisBufCoord(int x, int y, int ny) {
// return y*(ny+1) + x;
// }
/**
Transform cell-centered values into node-centered values
for visualization purposes
@param visBuffer pointer to visualization buffer
@param hd h-data on device
@param bd b-data on device
@param nx, ny x- and y-dimensions
@param dx, dy cellsize in x- and y-direction
@param average center (average) of h-data (water surface)
for visualization
@param scale scaling factor which is applied for visualization
@param offset offset factor for visualization
*/
__global__
void kernelCalcVisBuffer(float3* visBuffer, const float* hd, const float* bd,
int nx, int ny, float average, float scale, float offset)
{
int i = TILE_SIZE*blockIdx.x + threadIdx.x;
int j = TILE_SIZE*blockIdx.y + threadIdx.y;
if ((i <= nx) && (j <= ny)) {
int index = i*(nx+2) + j;
int index2 = (i+1)*(nx+2) + j;
visBuffer[j*(ny+1) + i] = make_float3(
i,
scaleFunction(0.25*(
hd[index]+hd[index + 1]+ hd[index2] + hd[index2 + 1] +
bd[index]+bd[index + 1]+ bd[index2] + bd[index2 + 1])
, average, scale, offset),
j);
}
//Corresponding C-Code:
// for (int j=0; j<ny+1;j++)
// for (int i=0;i<nx+1;i++)
// Vtk_file << i*dx<<" "<<j*dy<<" "
// << 0.25*(h[i][j]+h[i+1][j]+h[i][j+1]+h[i+1][j+1]
// +b[i][j]+b[i+1][j]+b[i][j+1]+b[i+1][j+1])
}
/**
copy h-values into visualization buffer
@param _visBuffer visualization buffer
*/
void Simulation::updateVisBuffer(float3* _visBuffer) {
const float* hd = splash->getCUDA_waterHeight();
const float* bd = splash->getCUDA_bathymetry();
int nx = splash->getNx();
int ny = splash->getNy();
// // Fill ghost layer corner cells
// kernelHdBufferEdges<<<1,1>>>(hd, nx, ny);
// Interpolate cell centered h-values
dim3 dimBlock(TILE_SIZE,TILE_SIZE);
dim3 dimGrid((nx+TILE_SIZE)/TILE_SIZE,(ny+TILE_SIZE)/TILE_SIZE);
kernelCalcVisBuffer<<<dimGrid,dimBlock>>>(_visBuffer, hd, bd, nx, ny,
wAverage, wScale, wOffset);
}
/**
Function used for debugging. Outputs the current visBuffer
@param visBuffer current visualization buffer
*/
void Simulation::debugVisBuffer(float3* _visBuffer) {
int nx = splash->getNx();
int ny = splash->getNy();
float* dbgBuf = new float[(nx+1)*(ny+1)*3];
int size = (nx+1)*(ny+1)*sizeof(float3);
cudaMemcpy(dbgBuf,_visBuffer, size, cudaMemcpyDeviceToHost);
checkCUDAError("memory of visualization buffer not transferred");
for (int i = 0; i < 3*(nx+1)*(ny+1); i = i + 3) {
std::cout << "(" << dbgBuf[i] << ", " << dbgBuf[i+1] << ", " << dbgBuf[i+2] << ") \n";
}
delete[] dbgBuf;
}
/**
Compute normal of a triangle
@param fVert1-fVert3 vertices of the triangle
@return resulting normal
*/
inline __device__
float3 calculateNormal(float3 fVert1, float3 fVert2,
float3 fVert3)
{
float Qx, Qy, Qz, Px, Py, Pz;
Qx = fVert2.x-fVert1.x;
Qy = fVert2.y-fVert1.y;
Qz = fVert2.z-fVert1.z;
Px = fVert3.x-fVert1.x;
Py = fVert3.y-fVert1.y;
Pz = fVert3.z-fVert1.z;
return make_float3(Py*Qz - Pz*Qy, Pz*Qx - Px*Qz, Px*Qy - Py*Qx);
}
/**
Compute new normals resulting from updated water surface
@param visBuffer vertex buffer holding water surface vertices
@param normalBuffer buffer which will contain new normals
@param nx, ny
*/
__global__
void kernelCalcNormals(float3* visBuffer, float3* normalBuffer, int nx, int ny)
{
int i = TILE_SIZE*blockIdx.x + threadIdx.x;
int j = TILE_SIZE*blockIdx.y + threadIdx.y;
if ((i < nx + 1) && (j < ny + 1)) {
int _i = i;
int _j = j;
// Handle boundaries
if (i >= nx) _i = i - 1;
if (j >= ny) _j = j - 1;
normalBuffer[j*(ny+1) + i] = calculateNormal(visBuffer[(_j*(ny+1) + _i)],
visBuffer[(_j+1)*(ny+1) + _i + 1],
visBuffer[(_j+1)*(ny+1) + _i]);
}
}
/**
Compute normals of the new water surface on graphics card
@param vertexBuffer cuda resource holding the vertex positions
@param destBuffer cuda resource holding the new normals
*/
void Simulation::calculateNormals(float3* vertexBuffer, float3* destBuffer) {
// splash->calculateNormals(vertexBuffer, destBuffer);
int nx = splash->getNx();
int ny = splash->getNy();
dim3 dimBlock(TILE_SIZE,TILE_SIZE);
dim3 dimGrid((nx+TILE_SIZE)/TILE_SIZE,(ny+TILE_SIZE)/TILE_SIZE);
kernelCalcNormals<<<dimGrid,dimBlock>>>(vertexBuffer, destBuffer, nx, ny);
}
/**
Helper function to check for errors in CUDA calls
param msg custom error message
source: NVIDIA
*/
#ifndef SIMULATION_H
#define SIMULATION_H
// =====================================================================
// This file is part of SWE_CUDA (see file SWE_Block.cu for details).
//
// Copyright (C) 2010,2011 Michael Bader, Kaveh Rahnema, Tobias Schnabel
//
// SWE_CUDA is free software: you can redristribute 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.
//
// SWE_CUDA 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 SWE_CUDA. If not, see <http://www.gnu.org/licenses/>.
// =====================================================================
#include <math.h>
#include <cuda_runtime.h>
#include "../SWE_BlockCUDA.hh"
#include "../scenarios/SWE_Scenario.h"
#include "../scenarios/SWE_VisInfo.h"
void checkCUDAError(const char *msg);
class Simulation {
public:
// Constructor + Destructor
Simulation (int nx, int ny, float dx, float dy,
SWE_Scenario* scene, SWE_VisInfo* visInfo, SWE_BlockCUDA* init_splash);
~Simulation();
// Restart simulation
void restart();
// Load new scenario after initialization
void loadNewScenario(SWE_Scenario* scene, SWE_VisInfo* visInfo);
// Save simulation state to file
void saveToFile();
// Return the bathymetry data
void setBathBuffer(float* output);
// Simulate single timestep on graphics card
void runCuda(struct cudaGraphicsResource **vbo_resource, struct cudaGraphicsResource **vbo_normals);
// Debugging
void writeDebugOutput(float3* destBuffer = NULL);
protected:
// Instance of SWE_BlockCUDA
SWE_BlockCUDA* splash;
// Current scenario
SWE_Scenario* myScenario;
// Current simulation time
float curTime;
// Is this our first simulation step?
int isFirstStep;
// Maximum of cell sizes
float maxCellSize;
// Scaling factors used by visualization
float bOffset, bScale, wOffset, wScale, bAverage, wAverage;
// Initalize boundaries defined by the scene
void initBoundaries(SWE_Scenario* scene);
// Initalize boundaries defined by an input file
int fileNumber;
// Use file input as boundary data
bool useFileInput;
// Holds maximum dimension
int maxDim;
// Compute new water surface
void calculateWaterSurface(float3* destBuffer);
// Compute normals of the water surface for shading
void calculateNormals(float3* vertexBuffer, float3* destBuffer);
void getScalingApproximation(const Float2D& h, const Float2D& b);
void updateVisBuffer(float3* _visBuffer);
void debugVisBuffer(float3* _visBuffer);
/**
Scale function used for visualization
*/
__device__ __host__
static float scaleFunction(float val, float average, float scale, float offset) {
return (val - average)*scale + offset;
}
static void calculateNormal(float fVert1[], float fVert2[],
float fVert3[], float fNormal[]);
};
#endif
// =====================================================================
// This file is part of SWE_CUDA (see file SWE_Block.cu for details).
//
// Copyright (C) 2010,2011 Tobias Schnabel
//
// SWE_CUDA is free software: you can redristribute 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.
//
// SWE_CUDA 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 SWE_CUDA. If not, see <http://www.gnu.org/licenses/>.
// =====================================================================
#include "visualization.h"
/**
Constructor. All dimensions are node-based, this means a grid consisting of
2x2 cells would have 3x3 nodes.
@param: window_title title of the window created
@param: _grid_x_size number of nodes of the grid (in x-direction)
@param: _grid_y_size number of nodes of the grid (in y-direction)
*/
Visualization::Visualization(int windowWidth, int windowHeight, const char* window_title,
int _grid_xsize, int _grid_ysize) {
// Initialize member variables
grid_xsize = _grid_xsize;
grid_ysize = _grid_ysize;
renderMode = SHADED;
// Initialize rendering
initSDL(windowWidth, windowHeight);
initGLWindow(windowWidth, windowHeight);
initGLDefaults();
initCUDA();
// Load camera and shaders
camera = new Camera(_grid_xsize*1.5f, window_title);
shaders = new Shader("vertex.glsl", "fragment.glsl");
if (shaders->shadersSupported()) {
printf("Shaders supported!\n");
} else {
printf("Shaders are NOT supported! Normal rendering mode\n");
}
if (shaders->shadersLoaded()) {
printf("Shaders successfully loaded!\n");
renderMode = WATERSHADER;
} else {
printf("Shaders error while loading shaders\n");
}
}
/**
Destructor (see note below)
*/
Visualization::~Visualization() {
delete camera;
delete shaders;
SDL_Quit();
}
/**
Frees all memory we used for geometry data
Needs to be called before destructor gets called
in order to work correctly
*/
void Visualization::cleanUp() {
deleteVBO(&vboWaterSurface, cuda_vbo_watersurface);
deleteVBO(&vboNormals, cuda_vbo_normals);
deleteVBO(&verticesIndex);
deleteVBO(&vboBathymetry);
}
/**
Main rendering function. Draws the scene and updates screen
*/
void Visualization::renderDisplay() {
// Set camera
camera->setCamera();
// Draw Scene
DrawBottom();
DrawBathymetry(vboBathymetry, verticesIndex);
// Shaded pass
DrawWaterSurface(vboWaterSurface, vboNormals, verticesIndex);
// Update framebuffer
camera->displayImage();
// Check for errors
if (glGetError() != GL_NO_ERROR) {
printf("OpenGL error occured..\n");
}
}
/**
Draws the water surface geometry (triangles)
@param vboID id of the vertex buffer object storing
the vertex positions
@param vboID id of the array buffer object storing
the corresponding normals for each vertex
@param verticesIndex id of the array buffer object storing the
vertex indices in rendering sequence
*/
void Visualization::DrawWaterSurface(GLuint vboID, GLuint vboNormals, GLuint verticesIndex) {
if (renderMode == WATERSHADER) {
// Depth pass first
glColorMask( GL_FALSE, GL_FALSE, GL_FALSE, GL_FALSE );
renderMode = SHADED;
DrawWaterSurface(vboID, vboNormals, verticesIndex);
renderMode = WATERSHADER;
glColorMask(GL_TRUE,GL_TRUE,GL_TRUE,GL_TRUE);
// Now render all visible geometry
shaders->enableShader();
glEnable(GL_BLEND);
glBlendFunc (GL_SRC_ALPHA, GL_ONE_MINUS_SRC_ALPHA);
}
else if (renderMode == SHADED) {
// Enable lighting
glEnable(GL_COLOR_MATERIAL);
glEnable(GL_LIGHTING);
} else {
glDisable(GL_LIGHTING);
}
// Enable array rendering
glEnableClientState(GL_NORMAL_ARRAY);
glEnableClientState(GL_VERTEX_ARRAY);
// Set rendering to VBO mode
glBindBuffer(GL_ARRAY_BUFFER,vboNormals);
glNormalPointer(GL_FLOAT, 0, 0);
glBindBuffer(GL_ARRAY_BUFFER, vboID);
glVertexPointer(3, GL_FLOAT, 0, 0);
glBindBuffer(GL_ELEMENT_ARRAY_BUFFER, verticesIndex);
// Enable VBO access and render triangles
glPushMatrix();
glTranslatef(-(grid_xsize-1)/2.0f,0.0f,-(grid_ysize-1)/2.0f);
glColor3f(0.3f*1.2f, 0.45f*1.2f, 0.9f*1.2f);
if (renderMode == WIREFRAME)
glPolygonMode( GL_FRONT_AND_BACK, GL_LINE );
glDrawElements(GL_TRIANGLES, 6*(grid_xsize - 1)*(grid_ysize - 1), GL_UNSIGNED_INT, NULL);
glPolygonMode( GL_FRONT_AND_BACK, GL_FILL );
glPopMatrix();
// Disable array rendering
glDisableClientState(GL_VERTEX_ARRAY);
glDisableClientState(GL_NORMAL_ARRAY);
glDisable(GL_LIGHTING);
glDisable(GL_BLEND);
shaders->disableShader();
}
/**
Draws the bathymetry geometry
@param vboID id of the vertex buffer object storing
the vertex positions
@param verticesIndex id of the array buffer object storing the
vertex indices in rendering sequence
*/
void Visualization::DrawBathymetry(GLuint vboID, GLuint verticesIndex) {
GLsizei stride = 6*sizeof(float);
// Enable lighting
glEnable(GL_COLOR_MATERIAL);
glEnable(GL_LIGHTING);
// Enable array rendering
glEnableClientState(GL_NORMAL_ARRAY);
glEnableClientState(GL_VERTEX_ARRAY);
// Bind buffers
glBindBuffer(GL_ARRAY_BUFFER, vboID);
glBindBuffer(GL_ELEMENT_ARRAY_BUFFER, verticesIndex);
// Set VBO pointers
const GLvoid* normalOffset = (GLvoid*) 12;
glNormalPointer(GL_FLOAT, stride, normalOffset);
glVertexPointer(3, GL_FLOAT, stride, 0);
// Render triangles
glPushMatrix();
glTranslatef(-(grid_xsize-1)/2.0f,0.0f,-(grid_ysize-1)/2.0f);
glColor3f(0.4f*1.1f, 0.36f*1.1f, 0.3f*1.1f);
glDrawElements(GL_TRIANGLES, 6*(grid_xsize - 1)*(grid_ysize - 1), GL_UNSIGNED_INT, NULL);
glPopMatrix();
// Disable array rendering
glDisableClientState(GL_VERTEX_ARRAY);
glDisableClientState(GL_NORMAL_ARRAY);
glDisable(GL_LIGHTING);
}
/**
Draws the bottom plane in our scene
*/
void Visualization::DrawBottom()
{
glBegin(GL_QUADS);
// Draw A Quad
glNormal3f(0.0f,1.0f,0.0f);
glColor3f(0.3f,0.3f,0.3f);
glVertex3f( (grid_xsize - 1)/2.0f,0.0f, (grid_ysize - 1)/2.0f); // Top Right Of The Quad (Bottom)
glVertex3f(-(grid_xsize - 1)/2.0f,0.0f, (grid_ysize - 1)/2.0f); // Top Left Of The Quad (Bottom)
glVertex3f(-(grid_xsize - 1)/2.0f,0.0f,-(grid_ysize - 1)/2.0f); // Bottom Left Of The Quad (Bottom)
glVertex3f( (grid_xsize - 1)/2.0f,0.0f,-(grid_ysize - 1)/2.0f); // Bottom Right Of The Quad (Bottom)
glEnd();
}
/**
Allocates memory for vertices and other geometry data.
@param sim instance of the simulation class
*/
void Visualization::init(Simulation* sim) {
// Create VBOs
createBathymetryVBO(&vboBathymetry, grid_xsize * grid_ysize * 6 * sizeof(float), sim);
createIndicesVBO(&verticesIndex, grid_xsize, grid_ysize);
createVertexVBO(&vboWaterSurface, grid_xsize * grid_ysize * 3 * sizeof(float),
&cuda_vbo_watersurface, cudaGraphicsMapFlagsNone);
createVertexVBO(&vboNormals ,grid_xsize * grid_ysize * 3 * sizeof(float),
&cuda_vbo_normals, cudaGraphicsMapFlagsWriteDiscard);
}
/**
Returns a pointer to the cuda memory object holding the
vertex normals
*/
cudaGraphicsResource** Visualization::getCudaNormalsPtr() {
return &cuda_vbo_normals;
}
/**
Returns a pointer to the cuda memory object holding the
vertex positions
*/
cudaGraphicsResource** Visualization::getCudaWaterSurfacePtr() {
return &cuda_vbo_watersurface;
}
/**
Returns, whether a special extension is supported by the current
graphics card
@param szTargetExtention string describing the extension to look for
*/
bool Visualization::IsExtensionSupported(const char* szTargetExtension )
{
const unsigned char *pszExtensions = NULL;
const unsigned char *pszStart;
unsigned char *pszWhere, *pszTerminator;
// Extension names should not have spaces
pszWhere = (unsigned char *) strchr( szTargetExtension, ' ' );
if( pszWhere || *szTargetExtension == '\0' )
return false;
// Get Extensions String
pszExtensions = glGetString( GL_EXTENSIONS );
// Search The Extensions String For An Exact Copy
pszStart = pszExtensions;
for(;;)
{
pszWhere = (unsigned char *) strstr( (const char *) pszStart, szTargetExtension );
if( !pszWhere )
break;
pszTerminator = pszWhere + strlen( szTargetExtension );
if( pszWhere == pszStart || *( pszWhere - 1 ) == ' ' )
if( *pszTerminator == ' ' || *pszTerminator == '\0' )
return true;
pszStart = pszTerminator;
}
return false;
}
/**
Initialize the simple directmedia layer (SDL) which is a wrapper library
for OpenGL
@param windowWidth width in pixels of the window to create
@param windowHeight height in pixels
*/
void Visualization::initSDL(int windowWidth, int windowHeight) {
// Initialize SDL system
if ( SDL_Init(SDL_INIT_VIDEO) < 0 ) {
fprintf(stderr, "Unable to initialize SDL: %s\n", SDL_GetError());
exit(1);
}
// Enable double buffering and center window
SDL_GL_SetAttribute( SDL_GL_DOUBLEBUFFER, 1 );
SDL_putenv(const_cast<char *>("SDL_VIDEO_CENTERED=center"));
SDL_WM_SetCaption("Loading...", NULL);
// Initialize window with OpenGL support
if ( SDL_SetVideoMode(windowWidth, windowHeight, 0, SDL_OPENGL | SDL_RESIZABLE | SDL_ANYFORMAT) == NULL ) {
fprintf(stderr, "Unable to create OpenGL screen: %s\n", SDL_GetError());
SDL_Quit();
exit(2);
}
// Check OpenGL extension(s)
if (!IsExtensionSupported( "GL_ARB_vertex_buffer_object" )) {
printf("Vertex Buffer Objects Extension not supported! Exit..\n");
SDL_Quit();
exit(1);
}
// Load Vertex Buffer Extension
glGenBuffers = (PFNGLGENBUFFERSARBPROC) SDL_GL_GetProcAddress("glGenBuffersARB");
glBindBuffer = (PFNGLBINDBUFFERARBPROC) SDL_GL_GetProcAddress("glBindBufferARB");
glBufferData = (PFNGLBUFFERDATAARBPROC) SDL_GL_GetProcAddress("glBufferDataARB");
glDeleteBuffers = (PFNGLDELETEBUFFERSARBPROC) SDL_GL_GetProcAddress("glDeleteBuffersARB");
}
/**
Initializes OpenGL projection and viewport settings
@param width width in pixels of the window to create
@param height height in pixels
*/
void Visualization::initGLWindow(int width, int height) {
GLfloat ratio;
// Protect against a divide by zero
if ( height == 0 )
height = 1;
ratio = ( GLfloat )width / ( GLfloat )height;
// Setup our viewport.
glViewport( 0, 0, ( GLint )width, ( GLint )height );
// change to the projection matrix and set our viewing volume.
glMatrixMode( GL_PROJECTION );
glLoadIdentity( );
float Zoom = 0.00005f;
glFrustum(-Zoom * width, Zoom * width, -Zoom * height, Zoom * height, 0.1f, 10.0f*grid_ysize);
// Switch back to the modelview
glMatrixMode( GL_MODELVIEW );
}
/**
Gets called when window gets resized
@param newWidth new window width in pixels
@param newHeight height in pixels
*/
int Visualization::resizeWindow(int newWidth, int newHeight) {
if ( SDL_SetVideoMode( newWidth, newHeight, 0,
SDL_OPENGL | SDL_RESIZABLE | SDL_ANYFORMAT ) == NULL )
{
fprintf( stderr, "Could not get a surface after resize: %s\n", SDL_GetError( ) );
return 1;
} else {
initGLWindow ( newWidth, newHeight );
return 0;
}
initGLWindow(newWidth, newHeight);
}
/**
Initializes various OpenGL settings
*/
void Visualization::initGLDefaults() {
// Set background color to white
glClearColor(1.0f, 1.0f, 1.0f, 0.0f);
// Smooth (interpolated) shading
glShadeModel( GL_SMOOTH );
// Reset the depth buffer
glClearDepth( 1.0f );
// Enable depth testing
glEnable( GL_DEPTH_TEST );
// Comparisons in depth buffer via <=
glDepthFunc( GL_LEQUAL );
// Some perspective corrections
glHint( GL_PERSPECTIVE_CORRECTION_HINT, GL_NICEST );
// Ensure that all of the following transformations will only
// apply to our modelview matrix (objects + camera)
glMatrixMode(GL_MODELVIEW);
glLoadIdentity();
// Lighting
GLfloat LightAmbient[]= { 0.0f, 0.0f, 0.0f, 1.0f };
GLfloat LightDiffuse[]= { 1.0f, 1.0f, 1.0f, 1.0f };
GLfloat LightPosition[]= { 1.0f, 1.0f, 0.0f, 0.0f };
glLightfv(GL_LIGHT0, GL_POSITION,LightPosition); // Position The Light
glLightfv(GL_LIGHT0, GL_AMBIENT, LightAmbient); // Setup The Ambient Light
glLightfv(GL_LIGHT0, GL_DIFFUSE, LightDiffuse); // Setup The Diffuse Light
glEnable(GL_LIGHT0);
glEnable(GL_NORMALIZE);
glDisable(GL_LIGHTING);
}
/**
Calculate 1-D array offset for 2-D coordinates
@param x x-position
@param y y-position
@param width size of a line
*/
int Visualization::coord(int x, int y, int width) {
if (width == -1) width = grid_xsize;
return (y*width + x);
}
/**
Creates a vertex buffer object in OpenGL
@param size size in bytes to allocate
@return vboID pointer to an integer depicting the id of the
new vertex buffer object
*/
void Visualization::createVertexVBO(GLuint* vboID, int size)
{
// Create 1 buffer object
glGenBuffers(1, vboID);
glBindBuffer(GL_ARRAY_BUFFER, *vboID);
// Initialize buffer object
glBufferData(GL_ARRAY_BUFFER, size, NULL, GL_DYNAMIC_DRAW);
// Switch to default buffer
glBindBuffer(GL_ARRAY_BUFFER, 0);
}
/**
Creates a vertex buffer object in OpenGL and loads it with the bathymetry
data from simulation
@param size size in bytes to allocate
@param sim pointer to an instance of the simulation class
@return vboID pointer to an integer depicting the id of the
new vertex buffer object
*/
void Visualization::createBathymetryVBO(GLuint* vboID, int size, Simulation* sim) {
GLfloat* vBathy = new GLfloat[size/sizeof(GLfloat)];
// Create buffer object for vertex indices
glGenBuffers(1, vboID);
glBindBuffer(GL_ARRAY_BUFFER, *vboID);
// Initialize buffer object
sim->setBathBuffer(vBathy);
glBufferData(GL_ARRAY_BUFFER, size, vBathy, GL_STATIC_DRAW);
// Switch to default buffer
glBindBuffer(GL_ARRAY_BUFFER, 0);
delete[] vBathy;
}
/**
Updates vertex buffer object with new bathymetry
data from simulation
@param sim pointer to an instance of the simulation class
*/
void Visualization::updateBathymetryVBO(Simulation* sim) {
int size = grid_xsize * grid_ysize * 6 * sizeof(float);
GLfloat* vBathy = new GLfloat[size/sizeof(GLfloat)];
// Select buffer
glBindBuffer(GL_ARRAY_BUFFER, vboBathymetry);
// Update buffer object
sim->setBathBuffer(vBathy);
glBufferData(GL_ARRAY_BUFFER, size, vBathy, GL_STATIC_DRAW);
// Switch to default buffer
glBindBuffer(GL_ARRAY_BUFFER, 0);
delete[] vBathy;
}
/**
Creates a vertex buffer object in OpenGL and an associated CUDA resource
@param size size in bytes to allocate
@param vbo_res_flags cuda flags for memory access
@return vboID pointer to an integer depicting the id of the
new vertex buffer object
@return vbo_res cuda structure created by this function
*/
void Visualization::createVertexVBO(GLuint* vboID, int size, struct cudaGraphicsResource **vbo_res,
unsigned int vbo_res_flags)
{
createVertexVBO(vboID, size);
cudaGraphicsGLRegisterBuffer(vbo_res, *vboID, vbo_res_flags);
checkCUDAError("Couldn't register GL buffer");
}
/**
Create an array buffer object which holds a list of vertex indices.
Groups of 3 consecutive vertices corresponding to the indices will be
rendered to a triangle.
This array buffer object therefore describes how single grid points (nodes)
get transformed into a triangle mesh (tesselation).
@param vbo pointer to an integer depicting the id of a
vertex buffer object
@param xsize number of grid nodes (in x-direction)
@params ysize number of grid nodes (in y-direction)
@return void
*/
void Visualization::createIndicesVBO(GLuint* vboID, int xsize, int ysize)
{
// Create an array describing the vertex indices to be drawn
xsize = xsize;
ysize = ysize;
int noVertices = (xsize-1)*(ysize-1)*6;
if ((xsize < 1) || (ysize < 1)) {
noVertices = 0;
}
GLuint* vIndices = new GLuint[noVertices];
for(int y=0; y < ysize - 1; y++)
{
for(int x=0; x < xsize - 1; x++)
{
// Create tessellation of the grid
vIndices[coord(x,y, xsize - 1)*6] = coord(x,y);
vIndices[coord(x,y, xsize - 1)*6 + 1] = coord(x+1,y+1);
vIndices[coord(x,y, xsize - 1)*6 + 2] = coord(x,y+1);
vIndices[coord(x,y, xsize - 1)*6 + 3] = coord(x,y);
vIndices[coord(x,y, xsize - 1)*6 + 4] = coord(x+1,y);
vIndices[coord(x,y, xsize - 1)*6 + 5] = coord(x+1,y+1);
}
}
// Create buffer object for vertex indices
glGenBuffers(1, vboID);
glBindBuffer(GL_ELEMENT_ARRAY_BUFFER, *vboID);
// Initialize buffer object
glBufferData(GL_ELEMENT_ARRAY_BUFFER, sizeof(GLuint)*noVertices, vIndices, GL_STATIC_DRAW);
// Switch to default buffer
glBindBuffer(GL_ELEMENT_ARRAY_BUFFER, 0);
delete[] vIndices;
}
/**
Frees memory used by a vertex buffer object
@param vbo pointer to an integer depicting the id of a
vertex buffer object
*/
void Visualization::deleteVBO(GLuint* vbo)
{
if (vbo) {
glBindBuffer(1, *vbo);
glDeleteBuffers(1, vbo);
*vbo = 0;
}
}
/**
Frees memory used by a vertex buffer object and a CUDA resource
@param vbo pointer to an integer depicting the id of a
vertex buffer object
@param vbo_res pointer to a CUDA resource structure
*/
void Visualization::deleteVBO(GLuint* vbo, struct cudaGraphicsResource *vbo_res)
{
if (vbo) {
cudaGraphicsUnregisterResource(vbo_res);
deleteVBO(vbo);
}
}
/**
Initialize CUDA device
*/
void Visualization::initCUDA() {
int driver, runtime;
cudaDeviceProp prop;
int dev;
// Display CUDA Versions
cudaDriverGetVersion(&driver);
cudaRuntimeGetVersion(&runtime);
printf("CUDA Driver Version: %d, Runtime Version: %d\n", driver, runtime);
// Select GPU for CUDA and OpenGL
memset( &prop, 0, sizeof( cudaDeviceProp ) );
prop.major = 1;
prop.minor = 0;
cudaChooseDevice( &dev, &prop );
cudaGLSetGLDevice( dev ) ;
}
/**
Sets current rendering mode
@param mode rendering mode
*/
void Visualization::setRenderingMode(RenderMode mode) {
renderMode = mode;
}
/**
Switches between 3 different rendering modes:
- Shaded: Use OpenGL shading
- Wireframe: Only render edges of each triangle
- Watershader: Use custom GLSL shader for water surface
*/
void Visualization::toggleRenderingMode() {
switch( renderMode ) {
case SHADED:
renderMode = WIREFRAME;
break;
case WIREFRAME:
// Skip watershader if shaders not loaded
if (shaders->shadersLoaded()) {
renderMode = WATERSHADER;
} else {
renderMode = SHADED;
}
break;
case WATERSHADER:
renderMode = SHADED;
break;
}
}
#ifndef VISUALIZATION_H
#define VISUALIZATION_H
// =====================================================================
// This file is part of SWE_CUDA (see file SWE_Block.cu for details).
//
// Copyright (C) 2010,2011 Tobias Schnabel
//
// SWE_CUDA is free software: you can redristribute 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.
//
// SWE_CUDA 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 SWE_CUDA. If not, see <http://www.gnu.org/licenses/>.
// =====================================================================
#include <SDL.h>
#include <SDL_opengl.h>
#include <cuda_runtime.h>
#include <cuda_gl_interop.h>
#include "camera.h"
#include "simulation.h"
#include "shader.h"
void checkCUDAError(const char *msg);
typedef enum RenderMode {
SHADED, WIREFRAME, WATERSHADER
} RenderMode;
class Visualization {
public:
// Constructor and Destructor
Visualization(int windowWidth, int windowHeight, const char* window_title, int _grid_xsize, int _grid_ysize);
~Visualization();
void init(Simulation* sim);
void cleanUp();
Camera* camera;
// Access to CUDA VBO pointers
cudaGraphicsResource** getCudaNormalsPtr();
cudaGraphicsResource** getCudaWaterSurfacePtr();
// Main rendering function
void renderDisplay();
// Helper functions
void setRenderingMode(RenderMode mode);
void updateBathymetryVBO(Simulation* sim);
void toggleRenderingMode();
int resizeWindow(int newWidth, int newHeight);
private:
// Init helper functions
void initSDL(int windowWidth, int windowHeight);
void initGLWindow(int width, int height);
void initGLDefaults();
void initCUDA();
bool IsExtensionSupported(const char* szTargetExtension );
// Drawing functions
void DrawWaterSurface(GLuint vboID, GLuint vboNormals, GLuint verticesIndex);
void DrawBathymetry(GLuint vboID, GLuint verticesIndex);
void DrawBottom();
int grid_xsize;
int grid_ysize;
// Vertex Buffer objects
GLuint vboBathymetry;
GLuint verticesIndex;
GLuint vboWaterSurface;
GLuint vboNormals;
struct cudaGraphicsResource* cuda_vbo_watersurface;
struct cudaGraphicsResource* cuda_vbo_normals;
// VBO management functions
void createIndicesVBO(GLuint* vboID, int xsize, int ysize);
void createVertexVBO(GLuint* vboID, int size);
void createBathymetryVBO(GLuint* vboID, int size, Simulation* sim);
void createVertexVBO(GLuint* vboID, int size, struct cudaGraphicsResource **vbo_res,
unsigned int vbo_res_flags);
void deleteVBO(GLuint* vbo);
void deleteVBO(GLuint* vbo, struct cudaGraphicsResource *vbo_res);
// Rendering mode
RenderMode renderMode;
// Shader helper class
Shader* shaders;
// Helper function
int coord(int x, int y, int width = -1);
// VBO Extension Function Pointers
PFNGLGENBUFFERSARBPROC glGenBuffers; // VBO Name Generation Procedure
PFNGLBINDBUFFERARBPROC glBindBuffer; // VBO Bind Procedure
PFNGLBUFFERDATAARBPROC glBufferData; // VBO Data Loading Procedure
PFNGLDELETEBUFFERSARBPROC glDeleteBuffers; // VBO Deletion Procedure
};
#endif
#ifndef __SWE_VISINFO_H
#define __SWE_VISINFO_H
// =====================================================================
// This file is part of SWE_CUDA (see file SWE_Block.cu for details).
//
// Copyright (C) 2010,2011 Michael Bader, Kaveh Rahnema, Tobias Schnabel
//
// SWE_CUDA is free software: you can redristribute 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.
//
// SWE_CUDA 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 SWE_CUDA. If not, see <http://www.gnu.org/licenses/>.
// =====================================================================
/**
* SWE_VisInfo defines an interface that can be used for online
* visualisation of a shallow water simulation.
* In particular, it provides information required for proper
* scaling of the involved variables.
*/
class SWE_VisInfo {
public:
// Scaling factors for custom visualization
// virtual bool useCustomScaling() { return false; };
virtual float waterHeightAtRest() { return 10.0f; };
virtual float waterDistanceFromGround() { return 10.0f; };
virtual float waterVerticalScaling() { return 10.0f; };
virtual float bathyVerticalCenter() { return 0.0f; };
virtual float bathyDistanceFromGround() { return 0.0f; };
virtual float bathyVerticalScaling() { return 10.0f; };
};
#endif
#include "SWE_VtkScenario.h"
#include <fstream>
#include <algorithm>
#include <stdlib.h>
#include <iostream>
#include <string>
#include <vector>
using namespace std;
/**
* Protected Constructor for class SWE_VtkScenario:
* called by readVtkFile to create an object of correct size
* NOTE: the provided number of cells in x and y direction includes ghost cells
* @param nx cells in x direction
* @param ny cells in y direction
*/
SWE_VtkScenario::SWE_VtkScenario(int nx, int ny)
: h(nx,ny),u(nx,ny),v(nx,ny),b(nx,ny)
{
cout << "Create VtkScenario of size " << nx << 'x' << ny << endl;
// set illegal value for dx and dy
// -> to be overwritten by readVtkFile
dx = -1.0; dy = -1.0;
}
/**
* get water height at specified position (piecewise constant interpolation);
* @param x position x-coordinate
* @param y position y-coordinate
*/
float SWE_VtkScenario::getWaterHeight(float x, float y)
{
if (x<0 || x >1) std::cout << "ERROR in x coordinate: " << x << endl << flush;
if (y<0 || y >1) std::cout << "ERROR in y coordinate: " << y << endl << flush;
int i = int((x-bndPos[BND_LEFT])/dx);
int j = int((y-bndPos[BND_BOTTOM])/dy);
if (i<0 || i >h.getCols())
std::cout << "ERROR in i coordinate: " << i
<< " with x coordinate: " << x
<< " offset " << (x-bndPos[BND_LEFT]) << '/' << dx << endl << flush;
if (j<0 || j >h.getRows()) std::cout << "ERROR in j coordinate: " << j
<< " with y coordinate: " << y << endl << flush;
return h[i][j];
}
/**
* get velocity (x-direction) at specified position (piecewise constant interpolation);
* @param x position x-coordinate
* @param y position y-coordinate
*/
float SWE_VtkScenario::getVeloc_u(float x, float y)
{
int i = int((x-bndPos[BND_LEFT])/dx);
int j = int((y-bndPos[BND_BOTTOM])/dy);
return u[i][j];
}
/**
* get velocity (x-direction) at specified position (piecewise constant interpolation);
* @param x position x-coordinate
* @param y position y-coordinate
*/
float SWE_VtkScenario::getVeloc_v(float x, float y)
{
int i = int((x-bndPos[BND_LEFT])/dx);
int j = int((y-bndPos[BND_BOTTOM])/dy);
return v[i][j];
}
/**
* get bathymetry at specified position (piecewise constant interpolation);
* @param x position x-coordinate
* @param y position y-coordinate
*/
float SWE_VtkScenario::getBathymetry(float x, float y)
{
int i,j;
if (x<0) {
std::cout << "request left of boundary" << x << endl << flush;
i = 0;
} else if (x >1) {
std::cout << "request right of boundary" << x << endl << flush;
i = b.getCols()-1;
} else {
i = int((x-bndPos[BND_LEFT])/dx);
};
if (y<0) {
std::cout << "request below boundary" << y << endl << flush;
j = 0;
} else if (y>1) {
std::cout << "request above boundary" << y << endl << flush;
j = b.getRows()-1;
} else {
j = int((y-bndPos[BND_BOTTOM])/dy);
};
return b[i][j];
}
/**
* Read VTK file which was written by SWE_Block.
* Use finite state automaton to parse input data.
* @param filename input filename
*/
SWE_VtkScenario* SWE_VtkScenario::readVtkFile(char const * filename)
{
// Helper variables
string line;
int state = 0;
int lineCount = 0;
int totalLines = 0;
int offsetCount = 0;
int dimx, dimy;
SWE_VtkScenario* scene;
float* h = NULL;
float* u = NULL;
float* v = NULL;
float* b = NULL;
int i,j;
int lastElem = 0;
size_t found;
vector<string> values;
// Open file
ifstream file(filename, ios::in);
if (!file.good()) {
std::cout << "Couldn't open file: " << filename << "\n";
return NULL;
}
// Count lines
totalLines = std::count(std::istreambuf_iterator<char>(file),
std::istreambuf_iterator<char>(), '\n');
file.seekg (0, ios::beg);
// Read line by line
while (std::getline(file,line)) {
// Skip empty lines
lineCount++;
if (! line.length()) continue;
// Convert to lower case
std::transform(line.begin(), line.end(),line.begin(), ::tolower);
// Output progress
if (lineCount % 10000 == 0 || lineCount == totalLines) {
std::cout << "\r" << "Reading input file: " <<
(lineCount*100)/totalLines << "% completed. " << flush;
}
switch(state)
{
// Header part
case 0:
found = line.find("# vtk datafile");
if (found!=string::npos) {
state = 1;
}
break;
case 1:
state = 2;
break;
case 2:
if (line.find("ascii") == 0) {
state = 3;
}
break;
case 3:
if (line.find("dataset structured_grid") == 0) {
state = 4;
}
break;
case 4:
found = line.find("dimensions");
if (found!=string::npos) {
stringExplode(line," ", &values);
if (values.size() != 4) {
state = -1;
} else {
dimx = atoi(values[1].c_str()) - 1;
dimy = atoi(values[2].c_str()) - 1;
cout << "Read grid dimensions " << dimx << 'x' << dimy << endl;
scene = new SWE_VtkScenario(dimx,dimy);
h = scene->h.elemVector();
u = scene->u.elemVector();
v = scene->v.elemVector();
b = scene->b.elemVector();
lastElem = dimx*dimy-1;
state = 5;
}
}
break;
// Data part
case 5:
// read grid data stored in VTK file
found = line.find(" 0 ");
if (found != string::npos && scene->dx <= 0.0f) {
stringExplode(line," ", &values);
scene->dx = (float) atof(values[0].c_str());
// TODO:
// --> this will change, if VTK files can have other origins than (0,0)
scene->bndPos[BND_BOTTOM] = 0.0;
scene->bndPos[BND_TOP] = dimx*scene->dx;
}
if (line.substr(0,2) == "0 " && scene->dy <= 0.0f) {
stringExplode(line," ", &values);
scene->dy = (float) atof(values[1].c_str());
// TODO:
// --> this will change, if VTK files can have other origins than (0,0)
scene->bndPos[BND_LEFT] = 0.0;
scene->bndPos[BND_RIGHT] = dimy*scene->dy;
}
found = line.find("scalars h");
if (found != string::npos) {
state = 10;
}
break;
case 10:
found = line.find("lookup_table");
if (found!=string::npos) {
state = 11;
offsetCount = 0;
}
break;
case 11:
// Read h-values
i = offsetCount%dimx; j= offsetCount/dimx;
h[i*dimy+j] = (float) atof(line.c_str());
if (offsetCount < lastElem)
offsetCount++;
else
state = 12;
break;
case 12:
found = line.find("scalars u");
if (found!=string::npos) {
state = 20;
}
break;
case 20:
found = line.find("lookup_table");
if (found!=string::npos) {
state = 21;
offsetCount = 0;
}
break;
case 21:
// Read u-values
i = offsetCount%dimx; j= offsetCount/dimx;
u[i*dimy+j] = (float) atof(line.c_str());
if (offsetCount < lastElem)
offsetCount++;
else
state = 22;
break;
case 22:
found = line.find("scalars v");
if (found!=string::npos) {
state = 30;
}
break;
case 30:
found = line.find("lookup_table");
if (found!=string::npos) {
state = 31;
offsetCount = 0;
}
break;
case 31:
// Read v-values
i = offsetCount%dimx; j= offsetCount/dimx;
v[i*dimy+j] = (float) atof(line.c_str());
if (offsetCount < lastElem)
offsetCount++;
else
state = 32;
break;
case 32:
found = line.find("scalars b");
if (found!=string::npos) {
state = 40;
}
break;
case 40:
found = line.find("lookup_table");
if (found!=string::npos) {
state = 41;
offsetCount = 0;
}
break;
case 41:
// Read b-values
i = offsetCount%dimx; j= offsetCount/dimx;
i = i*dimy+j;
b[i] = (float) atof(line.c_str());
h[i] = h[i] - b[i];
if (offsetCount < lastElem)
offsetCount++;
else
state = 42;
break;
default:
break;
}
}
std::cout << std::endl;
// Close file
file.close();
// delete allocated scene object
if (state >= 5 && state < 42)
delete scene;
// return NULL pointer in case of error
if (state < 42) {
std::cout << "Parsing Error " << state << " -> Omitting input file." << std::endl;
return NULL;
}
cout << "Read input data into scenario: " << endl
<< (*scene) << endl;
scene->computeHeightAtRest();
// return pointer to scene, if everything went OK
return scene;
}
/**
* Computes an estimate of the water height at rest
* (considering water height h plus bathymetry b)
* -> average value of h+b (considering only cells with h>0).
* Value is stored in member variable heightAtRest
*/
void SWE_VtkScenario::computeHeightAtRest()
{
float havg = 0.0;
int cnt = 0;
for(int i=0;i<h.getRows();i++)
for(int j=0;j<h.getCols();j++)
if (h[i][j] > 0.0f) {
havg += h[i][j] + b[i][j];
cnt++;
};
heightAtRest = havg/cnt;
}
/**
Split a string by seperator into substrings and return them as
an array
@param str string to split
@param separator delimiter string
@return results array containing all substrings
*/
void SWE_VtkScenario::stringExplode(string str, string separator, vector<string>* results){
size_t found;
results->clear();
found = str.find_first_of(separator);
while(found != string::npos){
if(found > 0){
results->push_back(str.substr(0,found));
}
str = str.substr(found+1);
found = str.find_first_of(separator);
}
if(str.length() > 0){
results->push_back(str);
}
}
//==================================================================
// external class-related methods
//==================================================================
/**
* overload operator<< such that data can be written via cout <<
* -> needs to be declared as friend to be allowed to access private data
*/
ostream& operator<<(ostream& os, SWE_VtkScenario& scene) {
int nx = scene.h.getCols();
int ny = scene.h.getRows();
os << "Gitterzellen: " << nx << "x" << ny << endl;
cout << "Wellenhoehe:" << endl;
for(int i=0; i<nx; i++) {
for(int j=0; j<ny; j++) {
os << scene.h[i][j] << " ";
};
os << endl;
};
cout << "Geschwindigkeit in x-Richtung:" << endl;
for(int i=0; i<nx; i++) {
for(int j=0; j<ny; j++) {
os << scene.u[i][j] << " ";
};
os << endl;
};
cout << "Geschwindigkeit in y-Richtung:" << endl;
for(int i=0; i<nx; i++) {
for(int j=0; j<ny; j++) {
os << scene.v[i][j] << " ";
};
os << endl;
};
cout << "Bathymetry:" << endl;
for(int i=0; i<nx; i++) {
for(int j=0; j<ny; j++) {
os << scene.b[i][j] << " ";
};
os << endl;
};
os << flush;
return os;
}
#ifndef __SWE_VTK_SCENARIO_H
#define __SWE_VTK_SCENARIO_H
// =====================================================================
// This file is part of SWE_CUDA (see file SWE_Block.cu for details).
//
// Copyright (C) 2010,2011 Michael Bader, Kaveh Rahnema, Tobias Schnabel
//
// SWE_CUDA is free software: you can redristribute 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.
//
// SWE_CUDA 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 SWE_CUDA. If not, see <http://www.gnu.org/licenses/>.
// =====================================================================
#include <iostream>
#include <stdio.h>
#include <fstream>
#include <string>
#include <vector>
#include "../tools/help.hh"
#include "SWE_Scenario.h"
using namespace std;
/**
* Scenario "VTK File":
* the object will read the values stored in a VTK file
* (i.e., use such a visualisation file as checkpoint for restarting a simulation)
*
* Note: - Boundary Conditions are NOT stored in a VTK file;
* OUTFLOW conditions are assumed at all boundaries
* - similar, 0.1 is always returnedby function endSimulation()
*
* TODO: Option to read a VTK container file
* (for example: read by different MPI processes, such that each MPI process
* will represent one of the container's VTK files)
*/
class SWE_VtkScenario : public SWE_Scenario {
public:
// SWE_VtkScenario does not offer a public constructor
// -> use this method to create a SWE_VtkScenario object
static SWE_VtkScenario* readVtkFile(char const * filename);
virtual float getWaterHeight(float x, float y);
virtual float getVeloc_u(float x, float y);
virtual float getVeloc_v(float x, float y);
virtual float getBathymetry(float x, float y);
virtual float waterHeightAtRest() { return heightAtRest; };
virtual void setWaterHeightAtRest(float whar) { heightAtRest = whar; };
virtual BoundaryType getBoundaryType(BoundaryEdge edge) { return OUTFLOW; };
virtual float getBoundaryPos(BoundaryEdge edge) { return bndPos[edge]; };
protected:
SWE_VtkScenario(int nx, int ny);
Float2D h;
Float2D u;
Float2D v;
Float2D b;
float heightAtRest;
float bndPos[4];
float dx;
float dy;
void computeHeightAtRest();
static void stringExplode(string str, string separator, vector<string>* results);
friend ostream& operator<< (ostream& os, SWE_VtkScenario& scene);
};
ostream& operator<< (ostream& os, const SWE_VtkScenario& scene);
#endif
#include "SWE_VtkScenarioVisInfo.h"
#include <fstream>
#include <algorithm>
#include <stdlib.h>
#include <iostream>
#include <string>
#include <vector>
#include <math.h>
using namespace std;
/**
* Protected Constructor for class SWE_VtkScenarioVisInfo:
* called by readVtkFile to create an object of correct size;
* only calls constructor of SWE_VtkScenario
* NOTE: the provided number of cells in x and y direction includes ghost cells
* @param nx cells in x direction
* @param ny cells in y direction
*/
SWE_VtkScenarioVisInfo::SWE_VtkScenarioVisInfo(int nx, int ny)
: SWE_VtkScenario(nx,ny)
{
}
/**
* Read VTK file which was written by SWE_Block.
* Use finite state automaton to parse input data.
* @param filename input filename
*/
SWE_VtkScenarioVisInfo* SWE_VtkScenarioVisInfo::readVtkFile(char const * filename)
{
// Helper variables
string line;
int state = 0;
int lineCount = 0;
int totalLines = 0;
int offsetCount = 0;
int dimx, dimy;
SWE_VtkScenarioVisInfo* scene;
float* h = NULL;
float* u = NULL;
float* v = NULL;
float* b = NULL;
size_t found;
vector<string> values;
// Open file
ifstream file(filename, ios::in);
if (!file.good()) {
std::cout << "Couldn't open file: " << filename << "\n";
return NULL;
}
// Count lines
totalLines = std::count(std::istreambuf_iterator<char>(file),
std::istreambuf_iterator<char>(), '\n');
file.seekg (0, ios::beg);
// Read line by line
while (std::getline(file,line)) {
// Skip empty lines
lineCount++;
if (! line.length()) continue;
// Convert to lower case
offsetCount++;
std::transform(line.begin(), line.end(),line.begin(), ::tolower);
// Output progress
if (lineCount % 10000 == 0 || lineCount == totalLines) {
std::cout << "\r" << "Reading input file: " <<
(lineCount*100)/totalLines << "% completed. " << flush;
}
switch(state)
{
// Header part
case 0:
found = line.find("# vtk datafile");
if (found!=string::npos) {
state = 1;
}
break;
case 1:
state = 2;
break;
case 2:
if (line.find("ascii") == 0) {
state = 3;
}
break;
case 3:
if (line.find("dataset structured_grid") == 0) {
state = 4;
}
break;
case 4:
found = line.find("dimensions");
if (found!=string::npos) {
stringExplode(line," ", &values);
if (values.size() != 4) {
state = -1;
} else {
dimx = atoi(values[1].c_str());
dimy = atoi(values[2].c_str());
scene = new SWE_VtkScenarioVisInfo(dimx,dimy);
h = scene->h.elemVector();
u = scene->u.elemVector();
v = scene->v.elemVector();
b = scene->b.elemVector();
state = 5;
}
}
break;
// Data part
case 5:
// read grid data stored in VTK file
found = line.find(" 0 ");
if (found != string::npos && scene->dx <= 0.0f) {
stringExplode(line," ", &values);
scene->dx = (float) atof(values[0].c_str());
// TODO:
// --> this will change, if VTK files can have other origins than (0,0)
scene->bndPos[BND_LEFT] = 0.0;
scene->bndPos[BND_RIGHT] = dimx*scene->dx;
}
if (line.substr(0,2) == "0 " && scene->dy <= 0.0f) {
stringExplode(line," ", &values);
scene->dy = (float) atof(values[1].c_str());
// TODO:
// --> this will change, if VTK files can have other origins than (0,0)
scene->bndPos[BND_BOTTOM] = 0.0;
scene->bndPos[BND_TOP] = dimy*scene->dy;
}
found = line.find("scalars h");
if (found != string::npos) {
state = 10;
}
break;
case 10:
found = line.find("lookup_table");
if (found!=string::npos) {
state = 11;
offsetCount = 0;
}
break;
case 11:
// Read h-values
h[offsetCount - 1] = (float) atof(line.c_str());
if (offsetCount >= (dimx - 1)*(dimy - 1)) {
state = 12;
}
break;
case 12:
found = line.find("scalars u");
if (found!=string::npos) {
state = 20;
}
break;
case 20:
found = line.find("lookup_table");
if (found!=string::npos) {
state = 21;
offsetCount = 0;
}
break;
case 21:
// Read u-values
u[offsetCount - 1] = (float) atof(line.c_str());
if (offsetCount >= (dimx - 1)*(dimy - 1)) {
state = 22;
}
break;
case 22:
found = line.find("scalars v");
if (found!=string::npos) {
state = 30;
}
break;
case 30:
found = line.find("lookup_table");
if (found!=string::npos) {
state = 31;
offsetCount = 0;
}
break;
case 31:
// Read v-values
v[offsetCount - 1] = (float) atof(line.c_str());
if (offsetCount >= (dimx - 1)*(dimy - 1)) {
state = 32;
}
break;
case 32:
found = line.find("scalars b");
if (found!=string::npos) {
state = 40;
}
break;
case 40:
found = line.find("lookup_table");
if (found!=string::npos) {
state = 41;
offsetCount = 0;
}
break;
case 41:
// Read b-values
b[offsetCount - 1] = (float) atof(line.c_str());
h[offsetCount - 1] = h[offsetCount - 1] - b[offsetCount - 1];
if (offsetCount >= (dimx - 1)*(dimy - 1)) {
state = 42;
}
break;
default:
break;
}
}
std::cout << std::endl;
// Close file
file.close();
// delete allocated scene object
if (state >= 5 && state < 42)
delete scene;
// return NULL pointer in case of error
if (state < 42) {
std::cout << "Parsing Error. Omitting input file." << std::endl;
return NULL;
}
// std::cout << "Compute scaling Info for read file\n" << flush;
scene->computeHeightAtRest();
scene->computeScalingApproximation();
// return pointer to scene, if everything went OK
return scene;
}
/**
Computes a first approximation of the scaling values needed
for visualization.
Gets called after variables have been read from VTK file;
determines the average, mininimum and maximum values of the
bathymetry and water surface data.
Uses latter values to estimate the scaling factors.
*/
void SWE_VtkScenarioVisInfo::computeScalingApproximation() {
// Averages of B and (B + H)
float avgB, avgH;
// Minimum values
float minB, minH;
// Maximum values
float maxB, maxH;
int nx = h.getRows();
int ny = h.getCols();
int maxDim = (nx > ny) ? nx : ny;
// first, compute averages:
avgH = 0.0f;
avgB = 0.0f;
int cnt = 0;
for(int i=0; i<nx; i++)
for(int j=0; j<ny; j++)
if (h[i][j]> 0) {
// Update averages
avgH += (h[i][j] + b[i][j]);
avgB += b[i][j];
cnt++;
};
bAverage = avgB/cnt;
hAverage = avgH/cnt;
// then, conpute min and max values, using average values as initial guess
minB = bAverage;
maxB = bAverage;
minH = hAverage;
maxH = hAverage;
for(int i=0; i<nx; i++)
for(int j=0; j<ny; j++)
if (h[i][j]> 0) {
float bij = b[i][j];
float hbij = h[i][j] + bij;
// Update minima
minH = (hbij < minH) ? hbij : minH;
minB = (bij < minB) ? bij : minB;
// Update maxima
maxH = (hbij > maxH) ? hbij : maxH;
maxB = (bij > maxB) ? bij : maxB;
};
std::cout << "Computed min and max values for visualisation info: \n"
<< "h: max: " << maxH << ", min: " << minH << endl
<< "b: max: " << maxB << ", min: " << minB << endl
<< flush;
if (fabs(maxB - minB) > 0.0001f) {
bScale = (maxDim/20.0f)/(maxB - minB);
} else {
bScale = (maxDim/20.0f);
}
bOffset = bScale*(bAverage - minB) + maxDim/15.0f;
// alternative values
bScale = minB;
bOffset = bAverage;
if ((maxH - minH) < 0.0001f) {
hScale = 1.0f/(maxH- minH);
} else {
hScale = 1.0f;
}
hScale = maxDim*(hScale);
hOffset = bOffset + (maxDim/5.0f);
bScale = minH;
hOffset = hAverage;
std::cout << "Computed values for visualisation info: "
<< "average h: " << hAverage << "( bOffset:" << hOffset << ", scale: " << hScale << ");"
<< "average b: " << bAverage << "( bOffset:" << bOffset << ", scale: " << bScale << ");"
<< endl << flush;
}
#ifndef __SWE_VTK_SCENARIO_VISINFO_H
#define __SWE_VTK_SCENARIO_VISINFO_H
// =====================================================================
// This file is part of SWE_CUDA (see file SWE_Block.cu for details).
//
// Copyright (C) 2010,2011 Michael Bader, Kaveh Rahnema, Tobias Schnabel
//
// SWE_CUDA is free software: you can redristribute 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.
//
// SWE_CUDA 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 SWE_CUDA. If not, see <http://www.gnu.org/licenses/>.
// =====================================================================
#include "../tools/help.hh"
#include "SWE_VtkScenario.h"
#include "SWE_VisInfo.h"
/**
* Scenario "VTK File" incl. info for visualisation:
* extends Scenario SWE_VtkScenario
*/
class SWE_VtkScenarioVisInfo : public SWE_VtkScenario, public SWE_VisInfo {
public:
// SWE_VtkScenario does not offer a public constructor
// -> use this method to create a SWE_VtkScenario object
static SWE_VtkScenarioVisInfo* readVtkFile(char const * filename);
// Scaling factors for custom visualization
virtual float waterHeightAtRest() { return hAverage; };
virtual float waterDistanceFromGround() { return hOffset; };
virtual float waterVerticalScaling() { return hScale; };
virtual float bathyVerticalCenter() { return bAverage; };
virtual float bathyDistanceFromGround() { return bOffset; };
virtual float bathyVerticalScaling() { return bScale; };
protected:
SWE_VtkScenarioVisInfo(int nx, int ny);
void computeScalingApproximation();
float bAverage;
float hAverage;
float bScale;
float hScale;
float bOffset;
float hOffset;
};
#endif
......@@ -41,7 +41,7 @@ class SWE_RadialDamBreakScenario : public SWE_Scenario {
public:
float getWaterHeight(float x, float y) {
return ( sqrt( (x-0.5f)*(x-0.5f) + (y-0.5f)*(y-0.5f) ) < 0.1f ) ? 12.0f: 10.0f;
return ( sqrt( (x-0.5f)*(x-0.5f) + (y-0.5f)*(y-0.5f) ) < 0.1f ) ? 10.5f: 10.0f;
};
virtual BoundaryType getBoundaryType(BoundaryEdge edge) { return OUTFLOW; };
......@@ -58,7 +58,8 @@ class SWE_BathymetryDamBreakScenario : public SWE_Scenario {
float getBathymetry(float x, float y) {
return ( std::sqrt( (x-500.f)*(x-500.f) + (y-500.f)*(y-500.f) ) < 50.f ) ? -250.f: -260.f;
};
virtual float endSimulation() { return (float) 15; };
virtual float endSimulation() { return (float) 15; };
virtual BoundaryType getBoundaryType(BoundaryEdge edge) { return OUTFLOW; };
......
#ifndef __SWE_SIMPLE_SCENARIOS_VIS_H
#define __SWE_SIMPLE_SCENARIOS_VIS_H
// =====================================================================
// This file is part of SWE_CUDA (see file SWE_Block.cu for details).
//
// Copyright (C) 2010,2011 Michael Bader, Kaveh Rahnema, Tobias Schnabel
//
// SWE_CUDA is free software: you can redristribute 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.
//
// SWE_CUDA 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 SWE_CUDA. If not, see <http://www.gnu.org/licenses/>.
// =====================================================================
#include <math.h>
#include "SWE_simple_scenarios.h"
#include "SWE_VisInfo.h"
/**
* Scenario "Radial Dam Break":
* elevated water in the center of the domain
*/
class SWE_RadialDamBreakScenarioVisInfo
: public SWE_RadialDamBreakScenario,
public SWE_VisInfo {
virtual float waterVerticalScaling() { return 5.0f; };
/* use inherited default implementations */
};
/**
* Scenario "Splashing Pool with custom scaling":
* intial water surface has a fixed slope (diagonal to x,y)
* shows how to use custom scaling to enhance visualization
* results
*/
class SWE_SplashingPoolScenarioVisInfo
: public SWE_SplashingPoolScenario,
public SWE_VisInfo {
public:
float waterDistanceFromGround() { return 9.0f; };
float waterVerticalScaling() { return 5.0f; };
float bathyVerticalCenter() { return 0.0f; };
float bathyDistanceFromGround() { return 0.0f; };
float bathyVerticalScaling() { return 0.0f; };
virtual BoundaryType getBoundaryType(BoundaryEdge edge) { return OUTFLOW; };
};
#endif
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment