demo_howto.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import os\n",
    "import numpy as np \n",
    "import damask"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%config Completer.use_jedi = False"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# some pre-requisites, after `source env/DAMASK.sh`\n",
    "damask_root = %env DAMASK_ROOT"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 1. input files (geometery, material, load)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## geometry file"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# geometry file (3d periodic grid)\n",
    "size = np.ones(3)*1e-5\n",
    "cells = [8,8,8]\n",
    "N_grains = 4\n",
    "seeds = damask.seeds.from_random(size,N_grains,cells)\n",
    "grid = damask.GeomGrid.from_Voronoi_tessellation(cells,size,seeds)\n",
    "grid.save(f'Polycrystal_{N_grains}_{cells[0]}x{cells[1]}x{cells[2]}')\n",
    "grid"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## material file"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "# material file has three sections (__homogenization__, __phase__, __material__)\n",
    "mat = damask.YAML(homogenization={}, phase={}, material={})"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "hom_single = {\n",
    "    'SX': {                            # (1) name of homogenization, SX->single crystal\n",
    "        'N_constituents': 1,           # (2) components number\n",
    "        'mechanical': {'type': 'pass'} # (3) type of homogenization\n",
    "    }\n",
    "}\n",
    "\n",
    "mat.update(homogenization=hom_single) # add to material.yaml file"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "phase_alu = {\n",
    "    'Aluminium':{                                    # (1) name of your phase\n",
    "        'lattice': 'cF',                             # (2) lattice specific\n",
    "        'mechanical': {                              # (3) mechanical properties\n",
    "            'output': ['F', 'P', 'F_e', 'L_p', 'O'], # Your desired output quantities\n",
    "            'elastic': {},\n",
    "            'plastic': {}\n",
    "        }\n",
    "    }\n",
    "}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "alu_elast = damask.YAML.load(os.path.join(damask_root,'examples/config/phase/mechanical/elastic/Hooke_Al.yaml'))\n",
    "alu_plast = damask.YAML.load(os.path.join(damask_root,'examples/config/phase/mechanical/plastic/phenopowerlaw_Al.yaml'))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "phase_alu['Aluminium']['mechanical']['elastic']['type'] = alu_elast['type']\n",
    "phase_alu['Aluminium']['mechanical']['elastic']['C_11'] = alu_elast['C_11']\n",
    "phase_alu['Aluminium']['mechanical']['elastic']['C_12'] = alu_elast['C_12']\n",
    "phase_alu['Aluminium']['mechanical']['elastic']['C_44'] = alu_elast['C_44']\n",
    "phase_alu['Aluminium']['mechanical']['plastic'] = dict(alu_plast)\n",
    "\n",
    "mat.update(phase=phase_alu) # add to material.yaml file"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "# material section -> grain orientation, components and etc\n",
    "## generate a simple random texture and add it to material section\n",
    "rnd = damask.Rotation.from_random(N_grains)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "mat_config = damask.ConfigMaterial(**mat)\n",
    "mat_config.pop('material') # material_add() requires no 'material' key exsists"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "# each material is single crystal (homogenization) and alu (phase)\n",
    "mat_config = mat_config.material_add(O=rnd,phase='Aluminium',homogenization='SX') \n",
    "mat_config.save('material.yaml')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## load file"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "## handy snippet to account for exclusiveness of stress/strain BC\n",
    "def make_P(F,fill=0):\n",
    "    return [make_P(i,fill) if isinstance(i,list) else\\\n",
    "            fill if i == 'x' else 'x' for i in F]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# each load file contains (1) solver option, (2) several loadsteps \n",
    "load_config = damask.LoadcaseGrid(\n",
    "    solver={'mechanical':'spectral_basic'}, # (1) solver\n",
    "    loadstep=[]                             # (2) loadsteps                     \n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# each load step specifies (1) BC (stress or strain), (2) time discretization, (3) output frequency\n",
    "## load step 1 (tensile in x) ##\n",
    "dot_F = [[1.e-3, 0 , 0 ],\n",
    "         [   0 ,'x', 0 ],\n",
    "         [   0 , 0 ,'x']]\n",
    "\n",
    "loadstep = {\n",
    "    'boundary_conditions':{              # (1) BC\n",
    "        'mechanical':{                \n",
    "            'dot_F':dot_F,\n",
    "            'P':make_P(dot_F)\n",
    "        }\n",
    "    },\n",
    "    'discretization':{'t':10.,'N':40},   # (2) time discretization\n",
    "    'f_out':4                            # (3) output frequency\n",
    "}\n",
    "\n",
    "load_config['loadstep'].append(loadstep)\n",
    "\n",
    "## load step 2 (tensile in x) ##\n",
    "dot_F = [[1.e-3, 0 , 0 ],\n",
    "         [   0 ,'x', 0 ],\n",
    "         [   0 , 0 ,'x']]\n",
    "\n",
    "loadstep = {\n",
    "    'boundary_conditions':{              # (1) BC\n",
    "        'mechanical':{                \n",
    "            'dot_F':dot_F,\n",
    "            'P':make_P(dot_F)\n",
    "        }\n",
    "    },\n",
    "    'discretization':{'t':40.,'N':40},   # (2) time discretization\n",
    "    'f_out':4                            # (3) output frequency\n",
    "}\n",
    "\n",
    "load_config['loadstep'].append(loadstep)\n",
    "load_config.save('load.yaml')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "right now you have a similar input file named as \n",
    "- `Polycrystal_{N_grains}_{cells}.vti`\n",
    "- `material.yaml`\n",
    "- `load.yaml`\n",
    "\n",
    "in your folder now, which resemble the ones in `ref/`\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 2. run simulation with DAMASK-grid solver"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fn_vti = f'Polycrystal_{N_grains}_{cells[0]}x{cells[1]}x{cells[2]}.vti'\n",
    "fn_load = 'load.yaml'\n",
    "fn_material = 'material.yaml'"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "%%time\n",
    "run_log = damask.util.run(f'DAMASK_grid -g {fn_vti} -l {fn_load}  -m {fn_material}')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(''.join(run_log))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 3. process result and visualize"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## basic post-processing"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fn_hdf = '_'.join((fn_vti.split('.')[0], fn_load.split('.')[0], fn_material.split('.')[0]))+'.hdf5'\n",
    "result = damask.Result(fn_hdf)\n",
    "result"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# post-processing: (1) add variables, (2) view certain increments, (3) export for paraview\n",
    "## (1) can act only on certain increments! when you have too large data\n",
    "result.add_stress_Cauchy()                      # Cauchy stress\n",
    "result.add_strain()                             # ln(V)    \n",
    "result.add_equivalent_Mises('sigma')            # von Mises value maybe more representative\n",
    "result.add_equivalent_Mises('epsilon_V^0.0(F)') # von Mises like strain based on ln(V) "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "## (2) various ways to view\n",
    "result.view(times=result.times_in_range(5,35))\n",
    "result.view(increments=-1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "## (3) various export options\n",
    "result.export_VTK()\n",
    "result.export_simulation_setup(target_dir='sim_setup')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## customized analysis"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "import seaborn as sns"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "grains = [0,1,3]\n",
    "grid   = damask.GeomGrid.load(fn_vti) # get back my grid\n",
    "# loop over all increments and store P(1,1) per grain and avg(F(1,1))\n",
    "data = {g:pd.DataFrame() for g in grains}\n",
    "for inc in result.get(['F','P']).values():\n",
    "    P = inc['P']\n",
    "    F = inc['F']\n",
    "    for g in grains:\n",
    "        points = grid.material.flatten(order='F')==g\n",
    "        P_11 = P[points,0,0].flatten()\n",
    "        F_11 = np.broadcast_to(np.average(F[:,0,0]),P_11.shape)\n",
    "        x = pd.DataFrame({'F_11':F_11,'P_11':P_11})\n",
    "        data[g] = pd.concat((data[g],x),ignore_index=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "for g in grains:\n",
    "   plot = sns.lineplot(y='P_11',x='F_11',data=data[g])\n",
    "\n",
    "fig = plot.get_figure()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## some in-line visualization with pyvista"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "need to do the following environment setup  (use _virtual environment_ or _conda_ )\n",
    "- `pip install 'pyvista[jupyter]>=0.38.1'` \n",
    "- `pip install ipywidgets 'pyvista[all,trame]'`\n",
    "- `pip install trame`\n",
    "\n",
    "_caution! there might be extra lib incompatibility!|_"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pyvista as pv\n",
    "pv.set_jupyter_backend('trame')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "data_grid = pv.read(fn_vti)\n",
    "\n",
    "pl = pv.Plotter()\n",
    "pl.add_mesh(data_grid, scalars='material',edge_color='k',show_edges=True)\n",
    "pl.add_axes()\n",
    "pl.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# get the last incr\n",
    "rv = result.view(increments=-1)\n",
    "rv.export_VTK()\n",
    "last_inc = result.increments[-1]\n",
    "fn_rv = fn_hdf.replace('.hdf5',f'_inc{last_inc}.vti')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# plot sigma_vM of last incr\n",
    "data_last = pv.read(fn_rv)\n",
    "\n",
    "pl = pv.Plotter()\n",
    "pl.add_mesh(data_last, scalars='phase/mechanical/sigma_vM / Pa', cmap='viridis')\n",
    "pl.add_axes()\n",
    "pl.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.10.12"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}