xref: /libCEED/examples/python/tutorial-5-operator.ipynb (revision cdf95791513f7c35170bef3ba2e19f272fe04533)
1{
2 "cells": [
3  {
4   "cell_type": "markdown",
5   "metadata": {},
6   "source": [
7    "# libCEED for Python examples\n",
8    "\n",
9    "This is a tutorial to illustrate the main feautures of the Python interface for [libCEED](https://github.com/CEED/libCEED/), the low-level API library for efficient high-order discretization methods developed by the co-design [Center for Efficient Exascale Discretizations](https://ceed.exascaleproject.org/) (CEED) of the [Exascale Computing Project](https://www.exascaleproject.org/) (ECP).\n",
10    "\n",
11    "While libCEED's focus is on high-order finite/spectral element method implementations, the approach is mostly algebraic and thus applicable to other discretizations in factored form, as explained in the [user manual](https://libceed.readthedocs.io/)."
12   ]
13  },
14  {
15   "cell_type": "markdown",
16   "metadata": {},
17   "source": [
18    "## Setting up libCEED for Python\n",
19    "\n",
20    "Install libCEED for Python by running"
21   ]
22  },
23  {
24   "cell_type": "code",
25   "execution_count": null,
26   "metadata": {},
27   "outputs": [],
28   "source": [
29    "! python -m pip install libceed"
30   ]
31  },
32  {
33   "cell_type": "markdown",
34   "metadata": {},
35   "source": [
36    "## CeedOperator\n",
37    "\n",
38    "Here we show some basic examples to illustrate the `libceed.Operator` class. In libCEED, a `libceed.Operator` defines the finite/spectral element operator associated to a `libceed.QFunction` (see [the API documentation](https://libceed.readthedocs.io/en/latest/libCEEDapi.html#finite-element-operator-decomposition))."
39   ]
40  },
41  {
42   "cell_type": "markdown",
43   "metadata": {},
44   "source": [
45    "* In the following example, we create and apply a CeedOperator for the mass matrix in 1D. By applying this operator to a vector of 1's, we compute the length of this 1D domain, similar to Ex1-Volume in the [tutorial-6-shell tutorial](./tutorial-6-shell.ipynb)"
46   ]
47  },
48  {
49   "cell_type": "code",
50   "execution_count": null,
51   "metadata": {},
52   "outputs": [],
53   "source": [
54    "import libceed\n",
55    "import numpy as np\n",
56    "\n",
57    "ceed = libceed.Ceed()\n",
58    "\n",
59    "nelem = 15\n",
60    "p = 5\n",
61    "q = 8\n",
62    "nx = nelem + 1\n",
63    "nu = nelem*(p-1) + 1\n",
64    "\n",
65    "# Vectors\n",
66    "x = ceed.Vector(nx)\n",
67    "x_array = np.zeros(nx)\n",
68    "for i in range(nx):\n",
69    "  x_array[i] = i / (nx - 1.0)\n",
70    "x.set_array(x_array, cmode=libceed.USE_POINTER)\n",
71    "\n",
72    "qdata = ceed.Vector(nelem*q)\n",
73    "u = ceed.Vector(nu)\n",
74    "v = ceed.Vector(nu)\n",
75    "\n",
76    "# Restrictions\n",
77    "indx = np.zeros(nx*2, dtype=\"int32\")\n",
78    "for i in range(nx):\n",
79    "  indx[2*i+0] = i\n",
80    "  indx[2*i+1] = i+1\n",
81    "rx = ceed.ElemRestriction(nelem, 2, 1, 1, nx, indx, cmode=libceed.USE_POINTER)\n",
82    "\n",
83    "indu = np.zeros(nelem*p, dtype=\"int32\")\n",
84    "for i in range(nelem):\n",
85    "  for j in range(p):\n",
86    "    indu[p*i+j] = i*(p-1) + j\n",
87    "ru = ceed.ElemRestriction(nelem, p, 1, 1, nu, indu, cmode=libceed.USE_POINTER)\n",
88    "strides = np.array([1, q, q], dtype=\"int32\")\n",
89    "rui = ceed.StridedElemRestriction(nelem, q, 1, q*nelem, strides)\n",
90    "\n",
91    "# Bases\n",
92    "bx = ceed.BasisTensorH1Lagrange(1, 1, 2, q, libceed.GAUSS)\n",
93    "bu = ceed.BasisTensorH1Lagrange(1, 1, p, q, libceed.GAUSS)\n",
94    "\n",
95    "# QFunctions\n",
96    "qf_setup = ceed.QFunctionByName(\"Mass1DBuild\")\n",
97    "qf_mass = ceed.QFunctionByName(\"MassApply\")\n",
98    "\n",
99    "# Setup operator\n",
100    "op_setup = ceed.Operator(qf_setup)\n",
101    "op_setup.set_field(\"dx\", rx, bx, libceed.VECTOR_ACTIVE)\n",
102    "op_setup.set_field(\"weights\", libceed.ELEMRESTRICTION_NONE, bx,\n",
103    "                   libceed.VECTOR_NONE)\n",
104    "op_setup.set_field(\"qdata\", rui, libceed.BASIS_COLLOCATED,\n",
105    "                   libceed.VECTOR_ACTIVE)\n",
106    "op_setup.check()\n",
107    "print('Setup operator: ', op_setup)\n",
108    "\n",
109    "# Mass operator\n",
110    "op_mass = ceed.Operator(qf_mass)\n",
111    "op_mass.set_field(\"u\", ru, bu, libceed.VECTOR_ACTIVE)\n",
112    "op_mass.set_field(\"qdata\", rui, libceed.BASIS_COLLOCATED, qdata)\n",
113    "op_mass.set_field(\"v\", ru, bu, libceed.VECTOR_ACTIVE)\n",
114    "op_mass.check()\n",
115    "print('Mass operator: ', op_mass)\n",
116    "\n",
117    "# Setup\n",
118    "op_setup.apply(x, qdata)\n",
119    "\n",
120    "# Apply mass matrix\n",
121    "u.set_value(1)\n",
122    "op_mass.apply(u, v)\n",
123    "\n",
124    "# Check\n",
125    "with v.array_read() as v_array:\n",
126    "  print('The length of the domain is l = %4.2f'%np.sum(v_array))"
127   ]
128  }
129 ],
130 "metadata": {
131  "kernelspec": {
132   "display_name": "Python 3",
133   "language": "python",
134   "name": "python3"
135  },
136  "language_info": {
137   "codemirror_mode": {
138    "name": "ipython",
139    "version": 3
140   },
141   "file_extension": ".py",
142   "mimetype": "text/x-python",
143   "name": "python",
144   "nbconvert_exporter": "python",
145   "pygments_lexer": "ipython3",
146   "version": "3.8.5"
147  }
148 },
149 "nbformat": 4,
150 "nbformat_minor": 4
151}
152