|
474 | 474 | "\n", |
475 | 475 | "State-space models are fundamental tools in control theory and signal processing that allow us to analyze complex dynamical systems by breaking them down into first-order differential equations. They are particularly important for structural dynamics problems because they can capture both the internal states (like position and velocity) and external influences (like forces) in a unified mathematical framework. A typical **state-space model** formulation might look like this:\n", |
476 | 476 | "\n", |
477 | | - "$$\n", |
478 | | - "x[k+1] \\sim \\mathcal{N}(A x[k] + B p[k], Q),\n", |
479 | | - "$$\n", |
480 | | - "$$\n", |
481 | | - "y[k] \\sim \\mathcal{N}(G x[k] + J p[k], R),\n", |
482 | | - "$$\n", |
| 477 | + "$$x[k+1] \\sim \\mathcal{N}(A x[k] + B p[k], Q),$$\n", |
| 478 | + "\n", |
| 479 | + "$$y[k] \\sim \\mathcal{N}(G x[k] + J p[k], R),$$\n", |
483 | 480 | "\n", |
484 | 481 | "where:\n", |
485 | 482 | "- $x[k]$ represents the system states at time-step $k$\n", |
|
506 | 503 | "cell_type": "markdown", |
507 | 504 | "metadata": {}, |
508 | 505 | "source": [ |
509 | | - "In this example, the dynamics of a structural system are governed by its **mass** ($ M $), **stiffness** ($ K $), and **damping** ($ C $) matrices, leading to the equation of motion:\n", |
| 506 | + "In this example, the dynamics of a structural system are governed by its **mass** ($M$), **stiffness** ($K$), and **damping** ($C$) matrices, leading to the equation of motion:\n", |
510 | 507 | "\n", |
511 | | - "$$\n", |
512 | | - "M \\ddot{x}(t) + C \\dot{x}(t) + K x(t) = p(t),\n", |
513 | | - "$$\n", |
| 508 | + "$$M \\ddot{x}(t) + C \\dot{x}(t) + K x(t) = p(t),$$\n", |
514 | 509 | "\n", |
515 | | - "where $ x(t) $ represents the displacements at each degree of freedom, and $ (t) $ is the external force applied to the system.\n", |
| 510 | + "where $x(t)$ represents the displacements at each degree of freedom, and $p(t)$ is the external force applied to the system.\n", |
516 | 511 | "\n", |
517 | 512 | "This model captures the essential dynamics of a multi-story structure while remaining computationally manageable. The system matrices are defined as follows:\n", |
518 | 513 | "\n", |
519 | | - "- $ M $ is the diagonal mass matrix representing the lumped masses at each floor, \n", |
520 | | - "- $ K $ is the stiffness matrix representing inter-floor lateral stiffness, and \n", |
521 | | - "- $ C $ is the proportional damping matrix reflecting energy dissipation." |
| 514 | + "- $M$ is the diagonal mass matrix representing the lumped masses at each floor, \n", |
| 515 | + "- $K$ is the stiffness matrix representing inter-floor lateral stiffness, and \n", |
| 516 | + "- $C$ is the proportional damping matrix reflecting energy dissipation." |
522 | 517 | ] |
523 | 518 | }, |
524 | 519 | { |
|
617 | 612 | "\n", |
618 | 613 | "We convert the structural system into its **discrete-time state-space form** for numerical simulation. Starting from the equation of motion:\n", |
619 | 614 | "\n", |
620 | | - "$$\n", |
621 | | - "M \\ddot{x}(t) + C \\dot{x}(t) + K x(t) = F(t),\n", |
622 | | - "$$\n", |
| 615 | + "$$M \\ddot{x}(t) + C \\dot{x}(t) + K x(t) = F(t),$$\n", |
623 | 616 | "\n", |
624 | 617 | "we introduce the state variable:\n", |
625 | | - "$$\n", |
626 | | - "z(t) = \\begin{bmatrix} x(t) \\\\ \\dot{x}(t) \\end{bmatrix},\n", |
627 | | - "$$\n", |
| 618 | + "\n", |
| 619 | + "$$z(t) = \\begin{bmatrix} x(t) \\\\ \\dot{x}(t) \\end{bmatrix},$$\n", |
| 620 | + "\n", |
628 | 621 | "which allows us to express the system as:\n", |
629 | | - "$$\n", |
630 | | - "\\dot{z}(t) = A_{\\text{c}} z(t) + B_{\\text{c}} p(t),\n", |
631 | | - "$$\n", |
| 622 | + "\n", |
| 623 | + "$$\\dot{z}(t) = A_{\\text{c}} z(t) + B_{\\text{c}} p(t),$$\n", |
| 624 | + "\n", |
632 | 625 | "where:\n", |
633 | | - "- $ A_{\\text{c}} = \\begin{bmatrix} 0 & I \\\\ -(M^{-1} K) & -(M^{-1} C) \\end{bmatrix} $\n", |
634 | | - "- $ B_{\\text{c}} = \\begin{bmatrix} 0 \\\\ M^{-1} S_p \\end{bmatrix} $\n", |
635 | | - "- $ S_p $ is the input selection matrix that determines where the external forces $ p(t) $ are applied.\n", |
636 | | - "\n", |
637 | | - "To perform simulations, the system is discretized using a time step $ \\Delta t $ as:\n", |
638 | | - "$$\n", |
639 | | - "z[k+1] = A z[k] + B p[k],\n", |
640 | | - "$$\n", |
| 626 | + "- $A_{\\text{c}} = \\begin{bmatrix} 0 & I \\\\ -(M^{-1} K) & -(M^{-1} C) \\end{bmatrix}$\n", |
| 627 | + "- $B_{\\text{c}} = \\begin{bmatrix} 0 \\\\ M^{-1} S_p \\end{bmatrix}$\n", |
| 628 | + "- $S_p$ is the input selection matrix that determines where the external forces $p(t)$ are applied.\n", |
| 629 | + "\n", |
| 630 | + "To perform simulations, the system is discretized using a time step $\\Delta t$ as:\n", |
| 631 | + "\n", |
| 632 | + "$$z[k+1] = A z[k] + B p[k],$$\n", |
| 633 | + "\n", |
641 | 634 | "where:\n", |
642 | | - "- $ A = e^{A_{\\text{c}} \\Delta t} $ is the **state transition matrix**.\n", |
643 | | - "- $ B = (A - I) A_{\\text{c}}^{-1} B_{\\text{c}} $ is the **input matrix**, obtained by integrating the continuous-time system.\n", |
| 635 | + "- $A = e^{A_{\\text{c}} \\Delta t}$ is the **state transition matrix**.\n", |
| 636 | + "- $B = (A - I) A_{\\text{c}}^{-1} B_{\\text{c}}$ is the **input matrix**, obtained by integrating the continuous-time system.\n", |
644 | 637 | "\n", |
645 | | - "This state-space representation forms the basis for propagating the system states during simulation.\n" |
| 638 | + "This state-space representation forms the basis for propagating the system states during simulation." |
646 | 639 | ] |
647 | 640 | }, |
648 | 641 | { |
|
692 | 685 | "source": [ |
693 | 686 | "### Generating Input Forces\n", |
694 | 687 | "\n", |
695 | | - "External forces $ p[k] $ acting on the system are modeled as **Gaussian white noise**:\n", |
| 688 | + "External forces $p[k]$ acting on the system are modeled as **Gaussian white noise**:\n", |
| 689 | + "\n", |
| 690 | + "$$p[k] \\sim \\mathcal{N}(\\mu, \\sigma^2),$$\n", |
696 | 691 | "\n", |
697 | | - "$$\n", |
698 | | - "p[k] \\sim \\mathcal{N}(\\mu, \\sigma^2),\n", |
699 | | - "$$\n", |
700 | | - "where $ \\mu $ is the mean and $ \\sigma $ controls the intensity of the force.\n", |
| 692 | + "where $\\mu$ is the mean and $\\sigma$ controls the intensity of the force.\n", |
701 | 693 | "\n", |
702 | | - "In this example, the inputs are generated independently at each time step $ k $ and across input channels to simulate random excitations, such as wind or seismic forces.\n", |
| 694 | + "In this example, the inputs are generated independently at each time step $k$ and across input channels to simulate random excitations, such as wind or seismic forces.\n", |
703 | 695 | "\n" |
704 | 696 | ] |
705 | 697 | }, |
|
735 | 727 | "\n", |
736 | 728 | "System responses, such as accelerations, are often measured at specific locations using sensors. The measurements are simulated using the equation:\n", |
737 | 729 | "\n", |
738 | | - "$$\n", |
739 | | - "y[k] = G x[k] + J p[k] + v[k],\n", |
740 | | - "$$\n", |
| 730 | + "$$y[k] = G x[k] + J p[k] + v[k],$$\n", |
| 731 | + "\n", |
741 | 732 | "where:\n", |
742 | | - "- $ G $ maps the system states $ x[k] $ to measured outputs.\n", |
743 | | - "- $ J $ maps the input forces $ p[k] $ to the measurements.\n", |
744 | | - "- $ v[k] \\sim \\mathcal{N}(0, \\sigma_y^2 I) $ is Gaussian noise representing sensor inaccuracies.\n", |
| 733 | + "- $G$ maps the system states $x[k]$ to measured outputs.\n", |
| 734 | + "- $J$ maps the input forces $p[k]$ to the measurements.\n", |
| 735 | + "- $v[k] \\sim \\mathcal{N}(0, \\sigma_y^2 I)$ is Gaussian noise representing sensor inaccuracies.\n", |
745 | 736 | "\n", |
746 | | - "The noise variance $ \\sigma_y^2 $ is chosen as a fraction of the true system response variance for realism.\n", |
| 737 | + "The noise variance $\\sigma_y^2$ is chosen as a fraction of the true system response variance for realism.\n", |
747 | 738 | "\n", |
748 | 739 | "In this example, **accelerations** are measured at selected degrees of freedom (e.g., nodes 1 and 4).\n" |
749 | 740 | ] |
|
804 | 795 | "\n", |
805 | 796 | "The structural response under applied forces is governed by the state-space equations:\n", |
806 | 797 | "\n", |
807 | | - "$$\n", |
808 | | - "\\begin{aligned}\n", |
| 798 | + "$$\\begin{aligned}\n", |
809 | 799 | "x[k+1] & = A x[k] + B p[k], \\\\\n", |
810 | 800 | "y[k] & = G_{\\text{full}} x[k] + J_{\\text{full}} p[k],\n", |
811 | | - "\\end{aligned}\n", |
812 | | - "$$\n", |
813 | | - "where $ x[k] $ are the system states, $ p[k] $ are the input forces, and $ y[k] $ are the **full-field responses**, i.e., the response at every degree of freedom in our structure.\n", |
| 801 | + "\\end{aligned}$$\n", |
| 802 | + "\n", |
| 803 | + "where $x[k]$ are the system states, $p[k]$ are the input forces, and $y[k]$ are the **full-field responses**, i.e., the response at every degree of freedom in our structure.\n", |
814 | 804 | "\n", |
815 | 805 | "\n", |
816 | 806 | "The function below returns:\n", |
817 | | - "- **True States**: $ x_{\\text{real}} $, propagated using $ A $ and $ B $.\n", |
818 | | - "- **Full-Field Responses**: $ y_{\\text{real}} $, incorporating both states and inputs.\n", |
819 | | - "- **Input Forces**: $ p_{\\text{real}} $, generated as stochastic excitations.\n", |
820 | | - "- **Response Matrices**: $ G_{\\text{full}} $ (state-to-response) and $ J_{\\text{full}} $ (input-to-response).\n", |
| 807 | + "- **True States**: $x_{\\text{real}}$, propagated using $ A $ and $ B $.\n", |
| 808 | + "- **Full-Field Responses**: $y_{\\text{real}}$, incorporating both states and inputs.\n", |
| 809 | + "- **Input Forces**: $p_{\\text{real}}$, generated as stochastic excitations.\n", |
| 810 | + "- **Response Matrices**: $G_{\\text{full}}$ (state-to-response) and $J_{\\text{full}}$ (input-to-response).\n", |
821 | 811 | "\n", |
822 | 812 | "These outputs simulate the physical behavior of the system and serve as the basis for inference. We keep the matrices because they will be used later when analyzing our results.\n" |
823 | 813 | ] |
|
881 | 871 | "source": [ |
882 | 872 | "### Augmented State-Space Model\n", |
883 | 873 | "\n", |
884 | | - "In structural health monitoring, external input forces $ p[k] $ acting on a structure, such as environmental loads or unknown excitations, are often not directly measurable. To estimate both the system states $ x[k] $ and these unknown input forces, we **augment the state vector** as follows:\n", |
| 874 | + "In structural health monitoring, external input forces $p[k]$ acting on a structure, such as environmental loads or unknown excitations, are often not directly measurable. To estimate both the system states $x[k]$ and these unknown input forces, we **augment the state vector** as follows:\n", |
885 | 875 | "\n", |
886 | | - "$$\n", |
887 | | - "\\tilde{x}[k] = \n", |
| 876 | + "$$\\tilde{x}[k] = \n", |
888 | 877 | "\\begin{bmatrix}\n", |
889 | 878 | "x[k] \\\\\n", |
890 | 879 | "p[k]\n", |
891 | | - "\\end{bmatrix}.\n", |
892 | | - "$$\n", |
| 880 | + "\\end{bmatrix}.$$\n", |
893 | 881 | "\n", |
894 | 882 | "This approach allows us to simultaneously infer the internal system states (e.g., displacements and velocities) and the unknown inputs using available measurements.\n", |
895 | 883 | "\n", |
896 | | - "\n", |
897 | 884 | "The augmented system dynamics are then expressed as:\n", |
898 | 885 | "\n", |
899 | | - "$$\n", |
900 | | - "\\begin{aligned}\n", |
| 886 | + "$$\\begin{aligned}\n", |
901 | 887 | "\\tilde{x}[k+1] & = A_{\\text{aug}} \\tilde{x}[k] + w[k], \\\\\n", |
902 | 888 | "y[k] & = G_{\\text{aug}} \\tilde{x}[k] + v[k],\n", |
903 | | - "\\end{aligned}\n", |
904 | | - "$$\n", |
| 889 | + "\\end{aligned}$$\n", |
905 | 890 | "\n", |
906 | 891 | "where:\n", |
907 | | - "- $ A_{\\text{aug}} $: Augmented state transition matrix. \n", |
908 | | - "- $ G_{\\text{aug}} $: Augmented measurement matrix. \n", |
909 | | - "- $ Q_{\\text{akf}} $: Augmented process noise covariance, capturing uncertainties in both states and inputs. \n", |
910 | | - "- $ w[k] $, $ v[k] $: Process and measurement noise. \n", |
| 892 | + "- $A_{\\text{aug}}$: Augmented state transition matrix. \n", |
| 893 | + "- $G_{\\text{aug}}$: Augmented measurement matrix. \n", |
| 894 | + "- $Q_{\\text{akf}}$: Augmented process noise covariance, capturing uncertainties in both states and inputs. \n", |
| 895 | + "- $w[k]$, $v[k]$: Process and measurement noise. \n", |
911 | 896 | "\n", |
912 | 897 | "##### Full-Field vs. Measurement Space \n", |
913 | 898 | "\n", |
914 | 899 | "To avoid confusion, we define two augmented measurement matrices: \n", |
915 | | - "- $ G_{\\text{aug}} $: Projects the augmented state vector $ \\tilde{x}[k] $ to the observed **sensor measurements** (e.g., accelerations at specific nodes). \n", |
916 | | - "- $ G^* $: The **augmented full-field measurement matrix**, which projects the augmented state vector to the full-field **system response**. This includes all degrees of freedom (displacements, velocities, and accelerations). \n", |
| 900 | + "- $G_{\\text{aug}}$: Projects the augmented state vector $\\tilde{x}[k]$ to the observed **sensor measurements** (e.g., accelerations at specific nodes). \n", |
| 901 | + "- $G^*$: The **augmented full-field measurement matrix**, which projects the augmented state vector to the full-field **system response**. This includes all degrees of freedom (displacements, velocities, and accelerations). \n", |
917 | 902 | "\n", |
918 | 903 | "The distinction is critical:\n", |
919 | | - "- $ G_{\\text{aug}} $ is used directly in the smoother to estimate states and inputs from limited measurements. \n", |
920 | | - "- $ G^* $ is used later to reconstruct the full response field for visualization and validation.\n", |
| 904 | + "- $G_{\\text{aug}}$ is used directly in the smoother to estimate states and inputs from limited measurements. \n", |
| 905 | + "- $G^*$ is used later to reconstruct the full response field for visualization and validation.\n", |
921 | 906 | "\n", |
922 | | - "For clarity, we will refer to the **augmented full-field matrix** as $ G^* $ throughout the rest of this example, whereas, in the code, this will be the `G_aug_fullfield` object.\n", |
| 907 | + "For clarity, we will refer to the **augmented full-field matrix** as $G^*$ throughout the rest of this example, whereas, in the code, this will be the `G_aug_fullfield` object.\n", |
| 908 | + "\n", |
| 909 | + "#### Noise Covariances \n", |
923 | 910 | "\n", |
924 | | - "#### Noise Covariances \n", |
925 | 911 | "In this step, the process and measurement noise covariances are assumed to be **known** or **pre-calibrated**. For example:\n", |
926 | | - "- The input force uncertainty $ Q_p $ is set to reflect significant variability. \n", |
927 | | - "- State noise covariance $ Q_x $ is chosen to reflect uncertainty in the model. \n", |
| 912 | + "- The input force uncertainty $Q_p$ is set to reflect significant variability. \n", |
| 913 | + "- State noise covariance $Q_x$ is chosen to reflect uncertainty in the model. \n", |
928 | 914 | "\n", |
929 | | - "The augmented noise covariance matrix $ Q_{\\text{akf}} $ combines these quantities:\n", |
| 915 | + "The augmented noise covariance matrix $Q_{\\text{akf}}$ combines these quantities:\n", |
930 | 916 | "\n", |
931 | | - "$$\n", |
932 | | - "Q_{\\text{akf}} =\n", |
| 917 | + "$$Q_{\\text{akf}} =\n", |
933 | 918 | "\\begin{aligned}\n", |
934 | 919 | "\\begin{bmatrix}\n", |
935 | 920 | "Q_x & 0 \\\\\n", |
936 | 921 | "0 & Q_p\n", |
937 | 922 | "\\end{bmatrix}\n", |
938 | | - "\\end{aligned}.\n", |
939 | | - "$$\n", |
940 | | - "\n" |
| 923 | + "\\end{aligned}.$$\n" |
941 | 924 | ] |
942 | 925 | }, |
943 | 926 | { |
|
1077 | 1060 | "\n", |
1078 | 1061 | "Here, we define our **Augmented Kalman Filter (AKF)** smoother using RxInfer. This probabilistic model estimates the system states and unknown input forces based on the measurements.\n", |
1079 | 1062 | "- **State Prior**: We start with a prior belief about the initial state, `x0`. \n", |
1080 | | - "- **State Transition**: At each time step, the system state evolves based on the transition matrix $ A $ and process noise covariance $ Q $:\n", |
1081 | | - " $$\n", |
1082 | | - " x[k] \\sim \\mathcal{N}(A x[k-1], Q).\n", |
1083 | | - " $$\n", |
| 1063 | + "- **State Transition**: At each time step, the system state evolves based on the transition matrix $A$ and process noise covariance $Q$:\n", |
| 1064 | + " $$x[k] \\sim \\mathcal{N}(A x[k-1], Q).$$\n", |
1084 | 1065 | "- **Measurements**: The observations (sensor data) are modeled as noisy measurements of the states:\n", |
1085 | | - " $$\n", |
1086 | | - " y[k] \\sim \\mathcal{N}(G x[k], R),\n", |
1087 | | - " $$\n", |
1088 | | - " where $ G $ maps the states to the measurements, and $ R $ is the measurement noise covariance.\n" |
| 1066 | + " $$y[k] \\sim \\mathcal{N}(G x[k], R),$$\n", |
| 1067 | + " where $G$ maps the states to the measurements, and $R$ is the measurement noise covariance.\n" |
1089 | 1068 | ] |
1090 | 1069 | }, |
1091 | 1070 | { |
|
1123 | 1102 | "\n", |
1124 | 1103 | "2. **Set the Initial State**: \n", |
1125 | 1104 | " We start with a prior belief about the first state, assuming it's zero with some process noise: \n", |
1126 | | - " $$\n", |
1127 | | - " x_0 \\sim \\mathcal{N}(0, Q_{\\text{akf}}).\n", |
1128 | | - " $$\n", |
| 1105 | + " $$x_0 \\sim \\mathcal{N}(0, Q_{\\text{akf}}).$$\n", |
1129 | 1106 | "\n", |
1130 | 1107 | "3. **Run the Smoother**: \n", |
1131 | 1108 | " We define a helper function to keep things tidy. This function calls RxInfer’s `infer` method, which does the heavy lifting for us.\n", |
|
1219 | 1196 | "\n", |
1220 | 1197 | "While the smoother estimates the system states, we often care about physical quantities like accelerations or displacements across the entire structure.\n", |
1221 | 1198 | "\n", |
1222 | | - "Using the **augmented full-field matrix** $ G^* $, we compute:\n", |
| 1199 | + "Using the **augmented full-field matrix** $G^*$, we compute:\n", |
1223 | 1200 | "- **Response means** from state means: \n", |
1224 | | - " $$\n", |
1225 | | - " \\mu_y[i] = G^* \\mu_x[i].\n", |
1226 | | - " $$ \n", |
| 1201 | + " $$\\mu_y[i] = G^* \\mu_x[i].$$ \n", |
1227 | 1202 | "- **Response uncertainties** from state covariances: \n", |
1228 | | - " $$\n", |
1229 | | - " \\sigma_y[i] = \\sqrt{\\text{diag}(G^* \\Sigma_x[i] {G^*}^\\top)}.\n", |
1230 | | - " $$ \n", |
| 1203 | + " $$\\sigma_y[i] = \\sqrt{\\text{diag}(G^* \\Sigma_x[i] {G^*}^\\top)}.$$ \n", |
1231 | 1204 | "\n", |
1232 | 1205 | "This gives us both the expected **responses** and their **uncertainties** at each time step. \n", |
1233 | 1206 | "\n", |
|
1255 | 1228 | " x_marginals,\n", |
1256 | 1229 | " G_aug_fullfield,\n", |
1257 | 1230 | " N_data::Int\n", |
1258 | | - " )\n", |
| 1231 | + ")\n", |
1259 | 1232 | " \n", |
1260 | 1233 | " # preallocate the full field response\n", |
1261 | 1234 | " y_means = Vector{Vector{Float64}}(undef, N_data) # vector of vectors\n", |
|
0 commit comments