<p style="font-family: Arial; font-size:3.75em;color:blue; font-style:bold">Dimensional Analysis</p>

## Find the relevant dimensionless parameters for a physics problem from a given set of physical variables in accordance with the Buckingham $\Pi$ theorem



```python
import numpy as np
import pandas as pd
import sympy as sp
from sympy import Matrix, init_printing

import itertools
```

#### The file ['physical-variables.tsv'](physical-variables.tsv) contains a table of physical variables and their asscociated dimensions (base quantity exponents)


```python
# Load physical variable dimensions table into pandas data frame
variables_df = pd.read_csv('physical-variables.tsv', delimiter='\t', index_col='Symbol')
variables_df.head()
```




<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }

    .dataframe tbody tr th {
        vertical-align: top;
    }

    .dataframe thead th {
        text-align: right;
    }
</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th>Quantity</th>
      <th>L</th>
      <th>M</th>
      <th>T</th>
      <th>I</th>
      <th>Theta</th>
      <th>N</th>
      <th>J</th>
      <th>SI unit</th>
      <th>in SI base units</th>
      <th>Category</th>
      <th>Comment</th>
    </tr>
    <tr>
      <th>Symbol</th>
      <th></th>
      <th></th>
      <th></th>
      <th></th>
      <th></th>
      <th></th>
      <th></th>
      <th></th>
      <th></th>
      <th></th>
      <th></th>
      <th></th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>l</th>
      <td>Length</td>
      <td>1</td>
      <td>0</td>
      <td>0</td>
      <td>0</td>
      <td>0</td>
      <td>0</td>
      <td>0</td>
      <td>m</td>
      <td>m</td>
      <td>Base quantity</td>
      <td>lowercase L</td>
    </tr>
    <tr>
      <th>m</th>
      <td>Mass</td>
      <td>0</td>
      <td>1</td>
      <td>0</td>
      <td>0</td>
      <td>0</td>
      <td>0</td>
      <td>0</td>
      <td>kg</td>
      <td>kg</td>
      <td>Base quantity</td>
      <td>-</td>
    </tr>
    <tr>
      <th>t</th>
      <td>Time</td>
      <td>0</td>
      <td>0</td>
      <td>1</td>
      <td>0</td>
      <td>0</td>
      <td>0</td>
      <td>0</td>
      <td>s</td>
      <td>s</td>
      <td>Base quantity</td>
      <td>-</td>
    </tr>
    <tr>
      <th>I</th>
      <td>Electric current</td>
      <td>0</td>
      <td>0</td>
      <td>0</td>
      <td>1</td>
      <td>0</td>
      <td>0</td>
      <td>0</td>
      <td>A</td>
      <td>A</td>
      <td>Base quantity</td>
      <td>-</td>
    </tr>
    <tr>
      <th>T</th>
      <td>Temperature</td>
      <td>0</td>
      <td>0</td>
      <td>0</td>
      <td>0</td>
      <td>1</td>
      <td>0</td>
      <td>0</td>
      <td>K</td>
      <td>K</td>
      <td>Base quantity</td>
      <td>-</td>
    </tr>
  </tbody>
</table>
</div>



#### Select the physical variables that are relevant to the physics problem! Choose variables (symbols) from the [physical-variables table](physical-variables.tsv)!


```python
# Enter the symbols separated by a space in a single string!

variables = 'F_D u eta rho D'              # Flow around sphere (drag force) 

#--------------------------------MORE EXAMPLES---------------------------------------
#variables = 'F_D eta u rho D'             # Flow around sphere (drag force)
#variables = 'v g rho t D eta m'           # velocity of a sphere sinking in a fluid
#variables = 'h u c_p rho eta k_th D'      # forced convection heat transfer
#variables = 'a_s p rho'                   # speed of sound
#variables = 'M_lambda k_B h_P c lambda T' # spectral exitance of black body
#variables = 'x F rho E_Y g L'             # large elastic deformations
#variables = 'F J E_Y L'                   # buckling load of a beam
#variables = 'f L E_Y rho  D'              # natural frequency of an elastic rod
#variables = 'M omega eta rho R'           # torque of rotating disc in fluid
#variables = 'P omega u L rho'             # power of wind turbine 
#variables = 'E p_0 R t rho_0'             # strong explosion
#-----------------------------------------------------------------------------------

# Split the string at the spaces to get a list of variables
variables_list  = variables.split()

# Display the selected variables and associated SI units
variables_table = variables_df.loc[variables_list,['Quantity', 
                                                   'SI unit',
                                                   'in SI base units',
                                                   'Category', 
                                                   'Comment']]
print('Selected physical variables: \n')
display(variables_table)

# Display the base quantity exponents for the selected variables
dimensions_table = variables_df.loc[variables_list,'L':'J'].transpose()
print('Dimensions: \n')
display(dimensions_table)
```

    Selected physical variables: 
    



<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }

    .dataframe tbody tr th {
        vertical-align: top;
    }

    .dataframe thead th {
        text-align: right;
    }
</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th>Quantity</th>
      <th>SI unit</th>
      <th>in SI base units</th>
      <th>Category</th>
      <th>Comment</th>
    </tr>
    <tr>
      <th>Symbol</th>
      <th></th>
      <th></th>
      <th></th>
      <th></th>
      <th></th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>F_D</th>
      <td>Drag / Resistance force</td>
      <td>N</td>
      <td>kg m / s^2</td>
      <td>Mechanics</td>
      <td>-</td>
    </tr>
    <tr>
      <th>u</th>
      <td>Flow velocity</td>
      <td>m / s</td>
      <td>m / s</td>
      <td>Kinematics</td>
      <td>-</td>
    </tr>
    <tr>
      <th>eta</th>
      <td>Dynamic viscosity</td>
      <td>Pa s</td>
      <td>kg / (m s)</td>
      <td>Mechanics</td>
      <td>-</td>
    </tr>
    <tr>
      <th>rho</th>
      <td>Density</td>
      <td>kg / m^3</td>
      <td>kg / m^3</td>
      <td>Mechanics</td>
      <td>-</td>
    </tr>
    <tr>
      <th>D</th>
      <td>Diameter</td>
      <td>m</td>
      <td>m</td>
      <td>Geometry</td>
      <td>-</td>
    </tr>
  </tbody>
</table>
</div>


    Dimensions: 
    



<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }

    .dataframe tbody tr th {
        vertical-align: top;
    }

    .dataframe thead th {
        text-align: right;
    }
</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th>Symbol</th>
      <th>F_D</th>
      <th>u</th>
      <th>eta</th>
      <th>rho</th>
      <th>D</th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>L</th>
      <td>1</td>
      <td>1</td>
      <td>-1</td>
      <td>-3</td>
      <td>1</td>
    </tr>
    <tr>
      <th>M</th>
      <td>1</td>
      <td>0</td>
      <td>1</td>
      <td>1</td>
      <td>0</td>
    </tr>
    <tr>
      <th>T</th>
      <td>-2</td>
      <td>-1</td>
      <td>-1</td>
      <td>0</td>
      <td>0</td>
    </tr>
    <tr>
      <th>I</th>
      <td>0</td>
      <td>0</td>
      <td>0</td>
      <td>0</td>
      <td>0</td>
    </tr>
    <tr>
      <th>Theta</th>
      <td>0</td>
      <td>0</td>
      <td>0</td>
      <td>0</td>
      <td>0</td>
    </tr>
    <tr>
      <th>N</th>
      <td>0</td>
      <td>0</td>
      <td>0</td>
      <td>0</td>
      <td>0</td>
    </tr>
    <tr>
      <th>J</th>
      <td>0</td>
      <td>0</td>
      <td>0</td>
      <td>0</td>
      <td>0</td>
    </tr>
  </tbody>
</table>
</div>


#### Base quantities (L, M, T, ... , J) with all zeros can be omitted


```python
# Drop all rows with all zeros
dimensions_table = dimensions_table[(dimensions_table.T != 0).any()]
display(dimensions_table)
table_index = list(dimensions_table.index)
```


<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }

    .dataframe tbody tr th {
        vertical-align: top;
    }

    .dataframe thead th {
        text-align: right;
    }
</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th>Symbol</th>
      <th>F_D</th>
      <th>u</th>
      <th>eta</th>
      <th>rho</th>
      <th>D</th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>L</th>
      <td>1</td>
      <td>1</td>
      <td>-1</td>
      <td>-3</td>
      <td>1</td>
    </tr>
    <tr>
      <th>M</th>
      <td>1</td>
      <td>0</td>
      <td>1</td>
      <td>1</td>
      <td>0</td>
    </tr>
    <tr>
      <th>T</th>
      <td>-2</td>
      <td>-1</td>
      <td>-1</td>
      <td>0</td>
      <td>0</td>
    </tr>
  </tbody>
</table>
</div>


#### The values of the dimensions table yield the dimensional matrix $A$

$
A=
\begin{pmatrix}
        a_{11}    & a_{12}  & ...     & a_{1n}\\
        a_{21}    & a_{22}  & ...     & a_{2n} \\
        \vdots    & \vdots  &         & \vdots   \\
        a_{m1}    & a_{m2}  & ...     & a_{mn}\\
\end{pmatrix}
$


```python
# Get dimensional matrix (values of dimensions_table)
A_mat = Matrix(dimensions_table.values)
print('A = ')
display(A_mat)

# Get the shape of the dimensional matrix
n_rows, n_cols = A_mat.shape
print(f'number of rows:    m = {n_rows}')
print(f'number of columns: n = {n_cols}')
```

    A = 



$\displaystyle \left[\begin{matrix}1 & 1 & -1 & -3 & 1\\1 & 0 & 1 & 1 & 0\\-2 & -1 & -1 & 0 & 0\end{matrix}\right]$


    number of rows:    m = 3
    number of columns: n = 5



#### The number of the expected dimensionless parameters can be obtained from the difference between the number of physical variables and the rank  of the dimensional matrix

#### $d=n-r$

$d$ : number of dimensionless parameters  
$n$ : number of physical variables  
$r$ : rank of dimensional matrix


```python
n_variables = n_cols
print(f'number of physical variables:       n = {n_variables}')

rank_A = A_mat.rank()
print(f'rank of dimensional matrix:         r = {rank_A}')

d = n_variables-rank_A
print(f'number of dimensionless parameters: d = {d}')

# Proceed only if there are solutions
if d > 0:
    proceed = True
else:
    proceed = False
    print('with the selected variables no dimensionless parameters can be found!')
```

    number of physical variables:       n = 5
    rank of dimensional matrix:         r = 3
    number of dimensionless parameters: d = 2


#### For this procedure to work, the $r$ rightmost column vectors of the dimensional matrix must be linearly independent. If this isn't the case, the variables of the input-string must be manually rearranged!


```python
if proceed:
    # Get submatrix C_temp containing the r rightmost vectors
    C_temp = Matrix(dimensions_table[variables_list[-rank_A:]].values)
    
    # Check if the rank of submatrix C_temp differs from rank_A
    if C_temp.rank() != rank_A:
        print(f'column-vectors of last {rank_A} variables are not linearly independent!')
        print(f'please rearrange the variables in the input-string!')
        print(f'possible combinations for the last {rank_A} variables are:')
        
        # Get all possible combinations of a number of rank_A variables
        combinations = itertools.combinations(variables_list, rank_A)
        
        # Check for each of the combinations if the rank of submatrix C_temp is equal to rank_A
        for combination in combinations:
            C_temp = Matrix(dimensions_table[list(combination)].values)
            if C_temp.rank() == rank_A:
                print('... ' + ' '.join(list(combination)))
                
        proceed = False                                                                   
```

#### If the number of columns of the dimensional matrix is lower or equal to the number of its rows ($n<= m$), or if the rank of the dimensional matrix is lower than the number of rows ($r<m$), then one or more rows must be deleted before we can proceed. However, the rank of the dimensional matrix must not change due to deletion of the rows!


```python
if proceed:
    if (n_cols <= n_rows) or (rank_A < n_rows):
        while (n_cols <= n_rows) or (rank_A < n_rows):
            print(' ')
            print('please delete a row!')
            
            # Calculate for each row the rank of the remaining matrix, if this row is deleted 
            ranks = []
            for i in range(n_rows):
                A_reduced = A_mat.copy()                               
                A_reduced.row_del(i)                                                                
                r_reduced = A_reduced.rank()                                                        
                ranks.append(r_reduced)                                                             
                
            # Create data frame with the ranks of the reduced matrices
            ranks_df = pd.DataFrame({'rank of remaining matrix, if this row is deleted':ranks,
                                     'this row may be deleted':(np.array(ranks)==rank_A),
                                     'row index':list(range(n_rows))},
                                     index=table_index)
            
            # Display dimensions table and ranks of reduced matrices in single data frame
            combined_dfs = dimensions_table.join(ranks_df)
            display(combined_dfs)
                        
            # Get user input with choice of row to be deleted
            row_to_del = np.NaN
            deletable_rows = list(ranks_df.loc[ranks_df['this row may be deleted'],'row index'])
            while row_to_del not in deletable_rows:                                                  
                print(f"please select suitable row index! Choose from {deletable_rows}")
                while True:
                    try:
                        row_to_del = int(input(f'Which row to delete? ... '))                       
                        break
                    except ValueError:
                        print("Oops!  That was no valid number.  Try again...")
                    
            # Delete row from matrix
            A_mat.row_del(row_to_del)
            n_rows, n_cols = A_mat.shape
            rank_A = A_mat.rank()
            print(rank_A)
            
            # Delete element from table index
            table_index.pop(row_to_del)
            
            # Delete row from dimension table
            dimensions_table.drop(dimensions_table.index[row_to_del], inplace=True)                  
            
        # Check if rightmost quadratic submatrix C_temp is invertible (i.e. not singular)
        C_temp = A_mat[:,-n_rows:]
        if C_temp.det() != 0:
            print('\n\nNew dimensions table')
            display(dimensions_table)
            print(f'number of rows: m = {n_rows}')
            print(f'number of columns: n = {n_cols}')
            
        else: 
            print(f'the rightmost quadratic submatrix is not invertible!')
            print(f'please rearrange the variables in the input-string!')
                       
            proceed = False 
```

#### Now that we have guaranteed that we have more columns than rows in the dimensional matrix ($n>m$), and that the matrix has full rank ($r=m$), we can start to solve following linear system of equations:

$
\begin{pmatrix}
        a_{11}    & a_{12}  & ...     & a_{1n}\\
        a_{21}    & a_{22}  & ...     & a_{2n} \\
        \vdots    & \vdots  &         & \vdots   \\
        a_{m1}    & a_{m2}  & ...     & a_{mn}\\
\end{pmatrix}
\cdot
\begin{pmatrix}
        k_{1}\\
        k_{2} \\
        \vdots  \\
        k_{n}\\
\end{pmatrix}
= 0
\hspace{1cm} \rightarrow \hspace{1cm} \textbf{A} \cdot \vec{k} = 0
$

where the $k_1, ..., k_n$ are the unknown exponentes we are looking for. Since we have more unknown variables ($n$) than equations ($m$), the linear system of equations is underdetermined. A homogeneous (with all constant terms equal to zero) underdetermined linear system always has non-trivial solutions. There are an infinity of such solutions, which form a vector space, whose dimension is the difference between the number of unknowns (in our case $n$) and the rank of the matrix of the system (in our case $r=m$). So our solutions vector space has the dimension $d=n-r$, and we can thus determine $d$ independent dimensionless parameters.

In order to find one of the possible solutions we can freely choose values for the first $d$ unknown variables ($k_1, ..., k_d$) and then solve the equation system for the remaining $r$ unknown variables ($k_{d+1}, ..., k_n$).

We therefore partition matrix $\textbf{A}$ into submatrix $\textbf{B}$ of size $(m \times d)$ and submatrix $\textbf{C}$ of size $(m \times (n-d))$. Since ($m=r$) and ($n-d=r$), $\textbf{C}$ is a quadratic matrix of size ($r \times r$). The $n$-elements vector $\vec{k}$ is split into the vector $\vec{x}$ with $d$ elements and vector $\vec{y}$ with $r$ elements. $\vec{x}$ contains the variables for which we want to set the values and $\vec{y}$ contains the variables of which we want to find the values. 

$
\begin{align}
\textbf{A} \cdot \vec{k} = 0 \\
\textbf{B} \cdot \vec{x} + \textbf{C} \cdot \vec{y} = 0
\end{align}
$

$
\begin{pmatrix}
        a_{11}    & a_{12}  & ...   & a_{1d}\\
        a_{21}    & a_{22}  & ...   & a_{2d} \\
         \vdots   & \vdots  &         & \vdots   \\
        a_{m1}    & a_{m2}  & ...   & a_{md}\\
\end{pmatrix}
\cdot
\begin{pmatrix}
        k_{1}\\
        k_{2} \\
        \vdots  \\
        k_{d}\\
\end{pmatrix}
+
\begin{pmatrix}
        a_{1(d+1)}    & a_{1(d+2)}  & ...   & a_{1n}\\
        a_{2(d+1)}    & a_{2(d+2)}  & ...   & a_{2n} \\
        \vdots        & \vdots      &       & \vdots   \\
        a_{m(d+1)}    & a_{m(d+2)}  & ...   & a_{mn}\\
\end{pmatrix}
\cdot
\begin{pmatrix}
        k_{d+1}\\
        k_{d+2} \\
        \vdots  \\
        k_{n}\\
\end{pmatrix}
= 0\\ 
$

We then subtract the first term to get it on the right side.

$
\textbf{C} \cdot \vec{y} = - \textbf{B} \cdot \vec{x}
$

$
\begin{pmatrix}
        a_{1(d+1)}    & a_{1(d+2)}  & ...   & a_{1n}\\
        a_{2(d+1)}    & a_{2(d+2)}  & ...   & a_{2n} \\
        \vdots        & \vdots      &       & \vdots   \\
        a_{m(d+1)}    & a_{m(d+2)}  & ...   & a_{mn}\\
\end{pmatrix}
\cdot
\begin{pmatrix}
        k_{d+1}\\
        k_{d+2} \\
        \vdots  \\
        k_{n}\\
\end{pmatrix}
= 
- \begin{pmatrix}
        a_{11}    & a_{12}  & ...   & a_{1d}\\
        a_{21}    & a_{22}  & ...   & a_{2d} \\
         \vdots   & \vdots  &         & \vdots   \\
        a_{m1}    & a_{m2}  & ...   & a_{md}\\
\end{pmatrix}
\cdot
\begin{pmatrix}
        k_{1}\\
        k_{2} \\
        \vdots  \\
        k_{d}\\
\end{pmatrix}\\
$

We then multiply on both sided of the equation from the left with the inverse of matrix $\textbf{C}$. On the left hand side the matrix multiplication of $\textbf{C}^{-1}$ with $\textbf{C}$ yields the identity matrix $\textbf{I}$ of size ($d \times d$) and the matrix multiplication $\textbf{I} \cdot \vec{y}$ yields of course the vector $\vec{y}$ itself.

$
\begin{align}
\textbf{C}^{-1} \cdot \textbf{C} \cdot \vec{y} &= - \textbf{C}^{-1} \cdot \textbf{B} \cdot \vec{x}\\
\textbf{I} \cdot \vec{y} &= - \textbf{C}^{-1} \cdot \textbf{B} \cdot \vec{x}
\end{align}
$

$
\begin{pmatrix}
        1      & 0      & ...   & 0\\
        0      & 1      & ...   & 0 \\
        \vdots & \vdots &       & \vdots   \\
        0      & 0      & ...   & 1\\
\end{pmatrix}
\cdot
\begin{pmatrix}
        k_{d+1}\\
        k_{d+2} \\
        \vdots  \\
        k_{n}\\
\end{pmatrix}
= 
- 
\begin{pmatrix}
        a_{1(d+1)}    & a_{1(d+2)}  & ...   & a_{1n}\\
        a_{2(d+1)}    & a_{2(d+2)}  & ...   & a_{2n} \\
        \vdots        & \vdots      &       & \vdots \\
        a_{m(d+1)}    & a_{m(d+2)}  & ...   & a_{mn}\\
\end{pmatrix}^{-1}
\cdot
\begin{pmatrix}
        a_{11}    & a_{12}  & ...   & a_{1d}\\
        a_{21}    & a_{22}  & ...   & a_{2d} \\
         \vdots   & \vdots  &         & \vdots   \\
        a_{m1}    & a_{m2}  & ...   & a_{md}\\
\end{pmatrix}
\cdot
\begin{pmatrix}
        k_{1}\\
        k_{2} \\
        \vdots  \\
        k_{d}\\
\end{pmatrix}\\
$

On the right hand side the term $- \textbf{C}^{-1} \cdot \textbf{B}$ yields the matrix $\textbf{D}$ of size $(m \times d)$.

$
\vec{y} = \textbf{D} \cdot \vec{x}
$

$
\begin{pmatrix}
        k_{d+1}\\
        k_{d+2} \\
        \vdots  \\
        k_{n}\\
\end{pmatrix}
= 
\textbf{D}
\cdot
\begin{pmatrix}
        k_{1}\\
        k_{2} \\
        \vdots  \\
        k_{d}\\
\end{pmatrix}\\
$

We can now freely choose the values for $\vec{x}$ and then determine the values of $\vec{y}$ to obtain all values needed to build the dimensionless parameters:

$\Pi = p_1^{k_{1}} \cdot p_2^{k_{2}} \cdot p_3^{k_{3}} \cdot \; ... \; \cdot p_n^{k_{n}}$

Since our solutions form a vector space with the dimenision of $d$, we can specify a different $\vec{x}$-vector $d$-times, but have to be sure, that these vectors are linearly independent. The easiest choice for the $\vec{x}$-vectors that meet this criteria are of course the unit vectors.

So, for the first dimensionless parameter we get:

$
\begin{pmatrix}
        k_{(1)d+1}\\
        k_{(1)d+2} \\
        \vdots  \\
        k_{(1)n}\\
\end{pmatrix}
= 
\textbf{D}
\cdot
\begin{pmatrix}
        1\\
        0 \\
        \vdots  \\
        0\\
\end{pmatrix}
\hspace{1cm} \rightarrow \hspace{1cm}
\Pi_1 = p_1^1 \cdot p_2^0 \cdot p_3^0 \cdot \; ... \cdot p_d^0 \cdot p_{d+1}^{k_{(1)d+1}} \cdot \; ... \; \cdot p_n^{k_{(1)n}}
$

for the second dimensionless parameter we get:

$
\begin{pmatrix}
        k_{(2)d+1}\\
        k_{(2)d+2} \\
        \vdots  \\
        k_{(2)n}\\
\end{pmatrix}
= 
\textbf{D}
\cdot
\begin{pmatrix}
        0\\
        1 \\
        \vdots  \\
        0\\
\end{pmatrix}
\hspace{1cm} \rightarrow \hspace{1cm}
\Pi_2 = p_1^0 \cdot p_2^1 \cdot p_3^0 \cdot \; ... \cdot p_d^0 \cdot p_{d+1}^{k_{(2)d+1}} \cdot \; ... \; \cdot p_n^{k_{(2)n}}
$

and so on...

However, we can simplify the procedure and obtain all values in one step. To do this, instead of $\vec{x}$ we enter all $d$ unit vectors at once (i.e., the identity matrix $\textbf{I}$ ) and obtain the solution matrix $\textbf{Y}$

$
\begin{align}
\vec{y}    &= \textbf{D} \cdot \vec{x}\\
\textbf{Y} &=\textbf{D} \cdot \textbf{I}\\
\textbf{Y} &=\textbf{D}\\
\end{align}
$

which shows that the matrix $\textbf{D}$ already contains all the values we are looking for!

To get a matrix $\textbf{K}$ with all the $k$-exponents we simply have to concatenate the identity matrix $\textbf{I}$ of size ($d \times d$) with the transpose of matrix $\textbf{D}$

$
\textbf{K} = \left[ I D^T \right]
$

$
\textbf{K} =
\begin{pmatrix}
        1       & 0      & ... & 0      & k_{(1)d+1} & k_{(1)d+2}  & ...   & k_{(1)n}\\
        0       & 1      & ... & 0      & k_{(2)d+1} & k_{(2)d+2}  & ...   & k_{(2)n} \\
        \vdots  & \vdots &     & \vdots & \vdots     & \vdots      &       & \vdots   \\
        0       & 0      & ... & 1      & k_{(d)d+1} & k_{(d)d+2}  & ...   & k_{(d)n}\\
\end{pmatrix}
$

where the i-th row of matrix $\textbf{K}$ contains the exponents of the i-th dimensionless parameter:

$
\Pi_i = p_1^{k_{(i)1}} \cdot p_2^{k_{(i)2}} \cdot p_3^{k_{(i)3}} \cdot \; ... \cdot p_n^{k_{(i)n}} \hspace{1cm} (i = 1 ...d)
$


```python
if proceed:
    print('A = ')
    display(A_mat)
    print('\n')
        
    # Partition matrix A into submatrix B...
    B_mat = A_mat[:,0:d]
    print('B = ')
    display(B_mat)
    print('\n')

    # ...and submatrix C
    C_mat = A_mat[:,d:]
    print('C = ')
    display(C_mat)
    print('\n')
    
    # Invert matrix C 
    C_inv = C_mat.inv()                         
    print('C_inv = ')
    display(C_inv)
    print('\n')
    
    # Calculate matrix D 
    D_mat = - C_inv*B_mat
    print('D = ')
    display(D_mat)
    print('\n')
    
    # Concatenate identity matrix with transpose of matrix D
    K_mat = sp.eye(d).row_join(D_mat.T)          
    print('K = ')
    display(K_mat)
    print('\n')
    
    # Create list of dimensionless parameter names
    Pi_list = []
    for j in range(d):
        Pi_list.append('Pi_' + str(j+1))
        
    # Create data frame from K_mat values
    exponents_df = pd.DataFrame(np.array(K_mat), index=Pi_list, columns=variables_list)
    display(exponents_df)
```

    A = 



$\displaystyle \left[\begin{matrix}1 & 1 & -1 & -3 & 1\\1 & 0 & 1 & 1 & 0\\-2 & -1 & -1 & 0 & 0\end{matrix}\right]$


    
    
    B = 



$\displaystyle \left[\begin{matrix}1 & 1\\1 & 0\\-2 & -1\end{matrix}\right]$


    
    
    C = 



$\displaystyle \left[\begin{matrix}-1 & -3 & 1\\1 & 1 & 0\\-1 & 0 & 0\end{matrix}\right]$


    
    
    C_inv = 



$\displaystyle \left[\begin{matrix}0 & 0 & -1\\0 & 1 & 1\\1 & 3 & 2\end{matrix}\right]$


    
    
    D = 



$\displaystyle \left[\begin{matrix}-2 & -1\\1 & 1\\0 & 1\end{matrix}\right]$


    
    
    K = 



$\displaystyle \left[\begin{matrix}1 & 0 & -2 & 1 & 0\\0 & 1 & -1 & 1 & 1\end{matrix}\right]$


    
    



<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }

    .dataframe tbody tr th {
        vertical-align: top;
    }

    .dataframe thead th {
        text-align: right;
    }
</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th>F_D</th>
      <th>u</th>
      <th>eta</th>
      <th>rho</th>
      <th>D</th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>Pi_1</th>
      <td>1</td>
      <td>0</td>
      <td>-2</td>
      <td>1</td>
      <td>0</td>
    </tr>
    <tr>
      <th>Pi_2</th>
      <td>0</td>
      <td>1</td>
      <td>-1</td>
      <td>1</td>
      <td>1</td>
    </tr>
  </tbody>
</table>
</div>


#### Now that we have our solution, we want to display the obtained dimensionless parameters nicely using Sympy.


```python
if proceed:
    # Define Sympy symbols
    var_symb = sp.symbols(variables)
    Pi_symb = sp.symbols(Pi_list)
    
    # Create list with Sympy equations
    expressions = []
    for i in range(d):
        product = 1
        for j in range(n_cols):
            product *= var_symb[j]**K_mat[i,j]
        expression = sp.Eq(Pi_symb[i],product)
        expressions.append(expression)
        display(expression)
        print('')

    display(variables_table)
```


$\displaystyle \Pi_{1} = \frac{F_{D} \rho}{\eta^{2}}$


    



$\displaystyle \Pi_{2} = \frac{D \rho u}{\eta}$


    



<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }

    .dataframe tbody tr th {
        vertical-align: top;
    }

    .dataframe thead th {
        text-align: right;
    }
</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th>Quantity</th>
      <th>SI unit</th>
      <th>in SI base units</th>
      <th>Category</th>
      <th>Comment</th>
    </tr>
    <tr>
      <th>Symbol</th>
      <th></th>
      <th></th>
      <th></th>
      <th></th>
      <th></th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>F_D</th>
      <td>Drag / Resistance force</td>
      <td>N</td>
      <td>kg m / s^2</td>
      <td>Mechanics</td>
      <td>-</td>
    </tr>
    <tr>
      <th>u</th>
      <td>Flow velocity</td>
      <td>m / s</td>
      <td>m / s</td>
      <td>Kinematics</td>
      <td>-</td>
    </tr>
    <tr>
      <th>eta</th>
      <td>Dynamic viscosity</td>
      <td>Pa s</td>
      <td>kg / (m s)</td>
      <td>Mechanics</td>
      <td>-</td>
    </tr>
    <tr>
      <th>rho</th>
      <td>Density</td>
      <td>kg / m^3</td>
      <td>kg / m^3</td>
      <td>Mechanics</td>
      <td>-</td>
    </tr>
    <tr>
      <th>D</th>
      <td>Diameter</td>
      <td>m</td>
      <td>m</td>
      <td>Geometry</td>
      <td>-</td>
    </tr>
  </tbody>
</table>
</div>


78f0f5be8b8c52ea4c7f90b02f958c9fc8df4208
ipynb
Jupyter Notebook
DimensionalAnalysis.ipynb
cwaidmann/Dimensional-Analysis
8b27df84137f45b934718aa2d448204e9cb4646e
MIT
DimensionalAnalysis.ipynb
cwaidmann/Dimensional-Analysis
8b27df84137f45b934718aa2d448204e9cb4646e
MIT
DimensionalAnalysis.ipynb
cwaidmann/Dimensional-Analysis
8b27df84137f45b934718aa2d448204e9cb4646e
MIT
Qwen/Qwen-72B
1. YES
2. YES