{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "31a94bbd-34ef-404b-a9be-4f3a216a8c3c",
   "metadata": {},
   "source": [
    "# Integration\n",
    "\n",
    "Here we show how to perform integration.\n",
    "\n",
    "This is a very simple example that simply used the attributes of the coef object"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f9ad6337-6ba4-47c7-b687-0bf119a5b637",
   "metadata": {},
   "source": [
    "#### Import general modules"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "8b4ca2be-46e8-443b-ab44-6d1a60abfec4",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Import required modules\n",
    "from mpi4py import MPI #equivalent to the use of MPI_init() in C\n",
    "import matplotlib.pyplot as plt\n",
    "import numpy as np\n",
    "\n",
    "# Get mpi info\n",
    "comm = MPI.COMM_WORLD"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "1702ba96-7854-4a2f-89d2-7f2b2aaa61b1",
   "metadata": {},
   "source": [
    "#### Import modules from pynek"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "09a4efe8-3557-4ff7-850e-0a26c8937bd6",
   "metadata": {},
   "outputs": [],
   "source": [
    "from pysemtools.io.ppymech.neksuite import pynekread\n",
    "from pysemtools.datatypes.msh import Mesh\n",
    "from pysemtools.datatypes.coef import Coef\n",
    "from pysemtools.datatypes.field import FieldRegistry"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "63764358",
   "metadata": {},
   "source": [
    "## Read the data and build objects\n",
    "\n",
    "In this instance, we create connectivity for the mesh object, given that we wish to use direct stiffness summation to reduce discontinuities."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "116a2e64",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "2024-08-25 21:02:37,603 - Mesh - INFO - Initializing empty Mesh object.\n",
      "2024-08-25 21:02:37,604 - Field - INFO - Initializing empty Field object\n",
      "2024-08-25 21:02:37,605 - pynekread - INFO - Reading file: ../data/rbc0.f00001\n",
      "2024-08-25 21:02:37,611 - Mesh - INFO - Initializing Mesh object from x,y,z ndarrays.\n",
      "2024-08-25 21:02:37,612 - Mesh - INFO - Initializing common attributes.\n",
      "2024-08-25 21:02:37,613 - Mesh - INFO - Creating connectivity\n",
      "2024-08-25 21:02:37,845 - Mesh - INFO - Mesh object initialized.\n",
      "2024-08-25 21:02:37,846 - Mesh - INFO - Mesh data is of type: float64\n",
      "2024-08-25 21:02:37,847 - Mesh - INFO - Elapsed time: 0.23524014099999999s\n",
      "2024-08-25 21:02:37,847 - pynekread - INFO - Reading field data\n",
      "2024-08-25 21:02:37,852 - pynekread - INFO - File read\n",
      "2024-08-25 21:02:37,853 - pynekread - INFO - Elapsed time: 0.248190697s\n",
      "2024-08-25 21:02:37,853 - Coef - INFO - Initializing Coef object\n",
      "2024-08-25 21:02:37,854 - Coef - INFO - Getting derivative matrices\n",
      "2024-08-25 21:02:37,856 - Coef - INFO - Calculating the components of the jacobian\n",
      "2024-08-25 21:02:37,919 - Coef - INFO - Calculating the jacobian determinant and inverse of the jacobian matrix\n",
      "2024-08-25 21:02:37,934 - Coef - INFO - Calculating the mass matrix\n",
      "2024-08-25 21:02:37,935 - Coef - INFO - Coef object initialized\n",
      "2024-08-25 21:02:37,936 - Coef - INFO - Coef data is of type: float64\n",
      "2024-08-25 21:02:37,936 - Coef - INFO - Elapsed time: 0.08297883800000005s\n"
     ]
    }
   ],
   "source": [
    "msh = Mesh(comm, create_connectivity=True)\n",
    "fld = FieldRegistry(comm)\n",
    "pynekread('../data/rbc0.f00001', comm, data_dtype=np.double, msh=msh, fld=fld)\n",
    "coef = Coef(msh, comm, get_area=False)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "45d0bf3a",
   "metadata": {},
   "source": [
    "## Integration in physical space\n",
    "\n",
    "In physical space, to integrate we must simply perform a weighted sum of a variable with the mass matrix. In this case, since we average, we divide by the mass matrix as well."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "022a0663",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Average w = -1.4410692531281081e-08\n"
     ]
    }
   ],
   "source": [
    "w_avrg = coef.glsum(fld.registry['w']*coef.B, comm, dtype=np.double)/coef.glsum(coef.B, comm, dtype=np.double)\n",
    "print('Average w =', w_avrg)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "23d815eb",
   "metadata": {},
   "source": [
    "The average of the z velocity in this case should be around 0"
   ]
  }
 ],
 "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.14"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}