/*
================================================================================
# BEGIN FILE: framer_fixed.c
# Date: June 27, 2025
# Author(s): Gary Argraves, Claude (Anthropic)
# 
# Purpose: Fixed C port with proper input parsing matching Pascal behavior
#
# Usage Example:
#   gcc -o framer_fixed framer_fixed.c -llapack -lblas -lm
#   ./framer_fixed < sample_1_story_3_bay.txt > output.html
#
# Key Fixes:
#   - Proper comment handling (ignore text after apostrophe)
#   - Correct input sequence matching Pascal version
#   - Skip blank lines appropriately
#   - Handle string input with embedded spaces
#   - Match exact Pascal input flow
#
================================================================================
*/

#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <math.h>
#include <time.h>
#include <ctype.h>
#include <stdbool.h>

/*=== LAPACK DECLARATIONS ===*/
extern void dgesv_(int* n, int* nrhs, double* a, int* lda, int* ipiv, 
                   double* b, int* ldb, int* info);

/*=== CONSTANTS ===*/
#define MAX_MATRIX_SIZE 50
#define MAX_MEMBERS 40
#define MAX_LOADS 50
#define MAX_MACROS 20
#define MAX_WIND_STORIES 10
#define MAX_MOMENTS 400
#define MAX_STRING 60
#define PUF 96

typedef char string60[MAX_STRING + 1];
typedef char string30[31];

/*=== GLOBAL VARIABLES (same as Pascal) ===*/
double d[MAX_MATRIX_SIZE + 1][MAX_MATRIX_SIZE + 1];
double b[MAX_MEMBERS + 1][7];
double c[MAX_MEMBERS + 1][3];
double e[MAX_WIND_STORIES + 1][5];
double ld[MAX_LOADS + 1][5];
int mp[MAX_MACROS + 1][3];
double r[MAX_MOMENTS + 1];

int s, p, q, z;
double l2, l3;
int lp, mc, lc;
int lm[MAX_LOADS + 1];
int psl, psm;

string60 proj_descrip, load_descrip, gWhat;
string30 time_str;
int option;
bool shear_flag;

int page_width = 1800, page_height = 1100;
int cur_x, cur_y, draw_pen, tab_cnt;
char plot_str[100000], tabulation_str[50000];
double gX, gY;
bool mouseover_enabled;

/*=== FIXED INPUT PARSING FUNCTIONS ===*/

void skip_comment_and_whitespace(void) {
    /* Skip whitespace and anything after apostrophe until newline */
    int c;
    while ((c = getchar()) != EOF) {
        if (c == '\'') {
            /* Skip rest of line after apostrophe */
            while ((c = getchar()) != EOF && c != '\n');
            if (c == '\n') break;
        } else if (c == '\n') {
            break;
        } else if (!isspace(c)) {
            /* Put back non-whitespace, non-comment character */
            ungetc(c, stdin);
            break;
        }
    }
}

double readr_fixed(void) {
    /* Read real number with proper comment and error handling */
    double r;
    int result;
    
    do {
        skip_comment_and_whitespace();
        result = scanf("%lf", &r);
        if (result != 1) {
            /* Skip to next line on error */
            int c;
            while ((c = getchar()) != EOF && c != '\n');
        }
    } while (result != 1);
    
    return r;
}

int readi_fixed(void) {
    /* Read integer with proper comment and error handling */
    int r;
    int result;
    
    do {
        skip_comment_and_whitespace();
        result = scanf("%d", &r);
        if (result != 1) {
            /* Skip to next line on error */
            int c;
            while ((c = getchar()) != EOF && c != '\n');
        }
    } while (result != 1);
    
    return r;
}

void reads_fixed(char *s, int max_len) {
    /* Read string until newline, handling comments properly */
    skip_comment_and_whitespace();
    
    if (fgets(s, max_len, stdin) != NULL) {
        /* Remove newline if present */
        char *newline = strchr(s, '\n');
        if (newline) *newline = '\0';
        
        /* Remove trailing whitespace */
        int len = strlen(s);
        while (len > 0 && isspace(s[len-1])) {
            s[--len] = '\0';
        }
    } else {
        s[0] = '\0';
    }
}

/*=== UTILITY FUNCTIONS ===*/
double pw(double n, int p) {
    double a = 1.0;
    for (int i = 1; i <= p; i++) a = a * n;
    return a;
}

void parse_up(char *s) {
    for (int i = 0; s[i]; i++) s[i] = toupper(s[i]);
}

void get_time_string(void) {
    time_t now;
    struct tm *tm_info;
    time(&now);
    tm_info = localtime(&now);
    strftime(time_str, sizeof(time_str), "%Y-%m-%d %H:%M:%S", tm_info);
}

void fatal_err(const char *s) {
    printf("%s\n", s);
    exit(1);
}

/*=== INITIALIZATION ===*/
void initial_var(void) {
    lp = 1;
    mc = 0;
    lc = 0;
    tab_cnt = 0;
    mouseover_enabled = false;
    
    for (int i = 0; i <= MAX_MATRIX_SIZE; i++) {
        for (int j = 0; j <= MAX_MATRIX_SIZE; j++) {
            d[i][j] = 0.0;
        }
    }
    
    for (int i = 0; i <= MAX_MEMBERS; i++) {
        for (int j = 0; j <= 6; j++) {
            b[i][j] = 0.0;
        }
        c[i][0] = c[i][1] = c[i][2] = 0.0;
    }
    
    for (int i = 0; i <= MAX_WIND_STORIES; i++) {
        for (int j = 0; j <= 4; j++) {
            e[i][j] = 0.0;
        }
    }
    
    for (int i = 0; i <= MAX_LOADS; i++) {
        for (int j = 0; j <= 4; j++) {
            ld[i][j] = 0.0;
        }
        lm[i] = 0;
    }
    
    for (int i = 0; i <= MAX_MACROS; i++) {
        mp[i][0] = mp[i][1] = mp[i][2] = 0;
    }
    
    strcpy(plot_str, "");
    strcpy(tabulation_str, "");
    strcpy(proj_descrip, "");
    strcpy(load_descrip, "");
    strcpy(gWhat, "");
}

void signon(void) {
    printf("RIGID FRAME ANALYSIS  By Gary Argraves - C Version 2025\n");
    printf("=======================================================\n");
    printf("Project Designation ? ");
    reads_fixed(proj_descrip, sizeof(proj_descrip));
    printf("%s\n", proj_descrip);
    get_time_string();
    printf("  Date: %s\n", time_str);
    printf("\n");
}

/*=== STRUCTURAL INPUT FUNCTIONS ===*/
void get_size(void) {
    printf("___INPUT FRAME SIZE___\n");
    printf("How many Stories ? ");
    s = readi_fixed();
    printf("%d\n", s);
    printf("How many Portals ? ");
    p = readi_fixed();
    printf("%d\n", p);
    printf("\n");
    q = (s + 1) * (p + 1) + s;
}

void get_geometry(void) {
    printf("___INPUT FRAME GEOMETRY (units=feet)___\n");
    
    /* Story heights input */
    printf("   Input story heights.\n");
    l3 = 0.0;
    for (int i = 1; i <= s; i++) {
        printf("Story #%2d ? ", i);
        double l = readr_fixed();
        printf("%8.2f\n", l);
        l3 += l;
        
        /* Assign height to all columns in this story */
        int z2 = (p + 1) * (s - i) + 1;
        for (int j = z2; j <= z2 + p; j++) {
            c[j][2] = l;
        }
    }
    printf("  Total height=%7.2f\n", l3);
    printf("\n");
    
    /* Bay widths input */
    printf("   Input beam lengths.\n");
    l2 = 0.0;
    for (int i = 1; i <= p; i++) {
        printf("Beam #%2d ? ", i);
        double l = readr_fixed();
        printf("%8.2f\n", l);
        l2 += l;
        
        /* Assign length to all beams in this bay across all stories */
        int j = i;
        do {
            b[j][2] = l;
            j += p;
        } while (j <= i + (s - 1) * p);
    }
    printf("  Total length=%7.2f\n", l2);
    
    /* Automatic scale calculation */
    psl = 5;
    psm = 50;
    while (l2 / psl > 9) {
        psl += 10;
        psm += 100;
    }

    /* Plotting parameters input */
    printf("\nPage width (pixels) =%5d ? ", page_width);
    page_width = readi_fixed();
    printf("\nPage height (pixels) =%5d ? ", page_height);
    page_height = readi_fixed();
    printf("\nPlot Scale Length (ft/in) =%3d ? ", psl);
    psl = readi_fixed();
    printf("\nPlot Scale Moment (kip-ft/inch) =%4d ? ", psm);
    psm = readi_fixed();
    printf("\n\n");
}

void get_member_sizes(void) {
    printf("___INPUT MEMBER SIZES - Moment of Inertia (I) (units=inch^4)___\n");
    
    for (int i = 1; i <= s; i++) {
        printf("\nSTORY #%2d\n", i);
        printf("-------\n");
        
        /* Beam sizes for this story */
        for (int j = 1; j <= p; j++) {
            printf("Beam #%2d ? ", j);
            b[j + p * (s - i)][1] = readr_fixed();
            printf("%10.2f\n", b[j + p * (s - i)][1]);
        }
        
        /* Column sizes for this story */
        for (int j = 1; j <= p + 1; j++) {
            printf("Column #%2d ? ", j);
            c[j + (s - i) * (p + 1)][1] = readr_fixed();
            printf("%10.2f\n", c[j + (s - i) * (p + 1)][1]);
        }
    }
    printf("\n");
}

void get_base_joint_type(void) {
    printf("___BASE JOINT TYPE___\n");
    printf("1. Fixed\n");
    printf("2. Pinned\n");
    do {
        printf("Enter 1 or 2 ? ");
        z = readi_fixed();
        printf("%d\n", z);
    } while (z != 1 && z != 2);
    printf("\n");
}

void get_macro_loads(void) {
    printf("  ___________ MACRO LOAD DEFINITION ____________\n");
    printf("\n");
    printf("            +-A--+---C---+\n");
    printf("         W  !    #########--W\n");
    printf("          \\ ====================== -- Story 2\n");
    printf("        -+-#I      I      I      I\n");
    printf("         | #I      I      I      I\n");
    printf("      C--| #I      I      I      I\n");
    printf("         | #====================== -- Story 1\n");
    printf("        -+-#I      I      I      I\n");
    printf("      A--|  I  P1  I  P2  I  P3  I\n");
    printf("       __|__I______I______I______I___ Ground\n");
    printf("\n");
    printf("  S=0 :Lateral Load   S=-1 :Next Macro Definition\n");
    printf("  S>0 :Gravity Load   S=-2 :Quit Macro Definition\n");
    printf("\n");

    int i;
    do {
        mc++;
        mp[mc][1] = lp;
        
        do {
            printf("Macro #%2d  Load #%d\n", mc, lp);
            printf("-------  ------\n");
            printf("Enter Story # S=? ");
            i = readi_fixed();
            printf("%d\n", i);
            
            if (i > s) {
                fatal_err("*** Error: Incorrect stories number ***");
            }
            
            double l4 = (i == 0) ? l3 : l2;
            
            if (i >= 0) {
                printf("Intensity W=? ");
                double w = readr_fixed();
                printf("%8.2f\n", w);
                
                printf("Location  A=? ");
                double j = readr_fixed();
                printf("%8.2f\n", j);
                
                printf("Length    C=? ");
                double k = readr_fixed();
                printf("%8.2f\n", k);
                
                if (j + k > l4) {
                    fatal_err("*** Error: A+C> Total Length of Frame ***");
                }
                
                ld[lp][1] = i;
                ld[lp][2] = w;
                ld[lp][3] = j;
                ld[lp][4] = k;
                mp[mc][2] = lp;
                lp++;
            } else if (mp[mc][2] == 0) {
                mc--;
            }
            printf("\n");
        } while (i >= 0);
    } while (i > -2);
}

/*=== LOAD PROCESSING (simplified version) ===*/
void process_lateral_load(double w, double j, double k) {
    double t = 0.0, t2 = 0.0;
    
    for (int h = 1; h <= s; h++) {
        int o = 1 + (s - h) * (p + 1);
        t = t + c[o][2];
        t2 = t - c[o][2];
        
        if ((j < t) && (j + k > t2)) {
            double a = 0.0;
            if (j > t2) a = j - t2;
            
            double b2 = c[o][2] - a;
            if (j + k < t) b2 = j + k - a - t2;
            
            double l = c[o][2];
            double c2 = a + b2;
            
            e[h][1] = e[h][1] - (w * fabs(1.0/12.0/pow(l,2) * 
                (6*pow(l,2)*(pow(c2,2)-pow(a,2)) - 8*l*(pw(c2,3)-pw(a,3)) + 3*(pw(c2,4)-pw(a,4)))));
            e[h][2] = e[h][2] + (w * fabs(1.0/12.0/pow(l,2) * 
                (4*l*(pw(c2,3)-pw(a,3)) - 3*(pw(c2,4)-pw(a,4)))));
            e[h][3] = e[h][3] + (w * b2 * (a + b2/2.0));
            e[h][4] = e[h][4] + (w * b2 * (l - a - b2/2.0));
        }
    }
}

void process_beam_load(int story, double w, double j, double k) {
    double t = 0.0, t2 = 0.0;
    
    for (int h = 1; h <= p; h++) {
        int o = h + (s - story) * p;
        t = t + b[o][2];
        t2 = t - b[o][2];
        
        if ((j < t) && (j + k > t2)) {
            double a = 0.0;
            if (j > t2) a = j - t2;
            
            double b2 = b[o][2] - a;
            if (j + k < t) b2 = j + k - a - t2;
            
            double l = b[o][2];
            double c2 = a + b2;
            
            b[o][3] = b[o][3] - (w * fabs(1.0/12.0/pow(l,2) * 
                (6*pow(l,2)*(pow(c2,2)-pow(a,2)) - 8*l*(pw(c2,3)-pw(a,3)) + 3*(pw(c2,4)-pow(a,4)))));
            b[o][4] = b[o][4] + (w * fabs(1.0/12.0/pow(l,2) * 
                (4*l*(pw(c2,3)-pw(a,3)) - 3*(pw(c2,4)-pw(a,4)))));
            b[o][5] = b[o][5] + (w * b2 * (a + b2/2.0));
            b[o][6] = b[o][6] + (w * b2 * (l - a - b2/2.0));
        }
    }
}

void apply_macro_loads(void) {
    printf("___APPLY LOADS TO FRAME___\n");
    printf("\n");
    printf("  # > 0 :Macro Loading\n");
    printf("  # < 0 :Single Load\n");
    printf("  # = 0 :Quit loading\n");
    printf("\n");
    
    lc = 0;
    printf("Enter Loading Designation ? ");
    reads_fixed(load_descrip, sizeof(load_descrip));
    parse_up(load_descrip);
    printf("%s\n", load_descrip);
    
    int m;
    do {
        printf("Enter # ? ");
        m = readi_fixed();
        printf("%d\n", m);
        
        if ((m > 0) && (m <= mc)) {
            for (int u = mp[m][1]; u <= mp[m][2]; u++) {
                lc++;
                lm[lc] = u;
                int i = (int)ld[u][1];
                double w = ld[u][2];
                double j = ld[u][3];
                double k = ld[u][4];
                
                if (i == 0) {
                    process_lateral_load(w, j, k);
                } else {
                    process_beam_load(i, w, j, k);
                }
            }
        } else if ((m < 0) && (m > -lp)) {
            int u = abs(m);
            lc++;
            lm[lc] = u;
            int i = (int)ld[u][1];
            double w = ld[u][2];
            double j = ld[u][3];
            double k = ld[u][4];
            
            if (i == 0) {
                process_lateral_load(w, j, k);
            } else {
                process_beam_load(i, w, j, k);
            }
        } else if (m != 0) {
            fatal_err("*** Error: No such Macro or Load. ***");
        }
    } while (m != 0);
    printf("\n");
}

/*=== MATRIX OPERATIONS (Full Pascal Implementation) ===*/
void setup_matrix(void) {
    /* Assemble global stiffness matrix using slope-deflection method - EXACT Pascal algorithm */
    printf("INFO: Setting up structural matrix...\n");
    
    int z2 = s + 1;
    if (z == 1) z2 = s;  /* Fixed base: reduce matrix size */
    
    /* Assembly follows exact Pascal algorithm */
    for (int u = 1; u <= z2; u++) {
        for (int j = s + u * (1 + p) - p; j <= s + u * (1 + p); j++) {
            
            /* Column below joint */
            if (u != 1) {
                double o = c[j - s - p - 1][1] / c[j - s - p - 1][2];
                if (j == s + u * (1 + p) - p) {
                    d[j][q + 1] = e[s + 2 - u][1];
                }
                d[j][j] = 4 * o + d[j][j];
                d[j][j - p - 1] = 2 * o;
                o = o / c[j - s - p - 1][2];
                d[j][u - 1] = 6 * o + d[j][u - 1];
                
                if (u == s + 1) goto NEXT_J;
                
                d[u][j] = -6 * o;
                d[u][j - p - 1] = -6 * o;
                d[u][u] = d[u][u] + 12 * o / c[j - s - p - 1][2];
                d[u][u - 1] = d[u][u - 1] - 12 * o / c[j - s - p - 1][2];
                d[j][u] = d[j][u] - 6 * o;
            }
            
            /* Column above joint */
            double o = c[j - s][1] / c[j - s][2];
            if (j == s + u * (1 + p) - p) {
                d[j][q + 1] = d[j][q + 1] + e[s + 1 - u][2];
            }
            d[j][j] = d[j][j] + 4 * o;
            
            if ((z != 1) || (u != s)) {
                d[u][j + p + 1] = 6 * o / c[j - s][2];
                d[j][j + p + 1] = 2 * o;
            }
            
            o = o / c[j - s][2];
            d[u][j] = d[u][j] + 6 * o;
            d[j][u] = d[j][u] + 6 * o;
            d[u][u] = d[u][u] + 12 * o / c[j - s][2];
            
            if (u != s) {
                d[j][u + 1] = d[j][u + 1] - 6 * o;
                d[u][u + 1] = d[u][u + 1] - 12 * o / c[j - s][2];
            }
            
            /* Beam to right of joint */
            if (j != s + u * (1 + p)) {
                o = b[j - s - u + 1][1] / b[j - s - u + 1][2];
                d[j][j + 1] = 2 * o;
                d[j][j] = d[j][j] + 4 * o;
                d[j][q + 1] = d[j][q + 1] + b[j - s - u + 1][3];
            }
            
            /* Beam to left of joint */
            if (j != s + u * (1 + p) - p) {
                o = b[j - s - u][1] / b[j - s - u][2];
                d[j][j] = d[j][j] + 4 * o;
                d[j][j - 1] = d[j][j - 1] + 2 * o;
                d[j][q + 1] = d[j][q + 1] + b[j - s - u][4];
            }
            
            NEXT_J:;
        }
    }

    /* Assemble displacement equilibrium equations */
    for (int i = 1; i <= s; i++) {
        int u = s + 1 - i;
        d[i][q + 1] = (e[u][3] + e[u][2] + e[u][1]) / c[(i - 1) * (p + 1) + 1][2];
        if (i > 1) {
            d[i][q + 1] = d[i][q + 1] + 
                (e[u + 1][4] - e[u + 1][2] - e[u + 1][1]) / c[(i - 2) * (p + 1) + 1][2];
        }
    }
    
    printf("INFO: Matrix assembly completed\n");
}

void solve_matrix(void) {
    /* Solve system using Pascal's Gaussian elimination algorithm */
    printf("INFO: Solving matrix system using Gaussian elimination...\n");
    
    bool matrix_err = false;
    bool restart;
    
    do {
        restart = false;
        int i = 1;
        
        do {
            /* Special handling for fixed base conditions */
            if ((z == 1) && (i == q - p)) break;
            
            /* Check for zero pivot */
            if (d[i][i] == 0.0) {
                if (i == q) {
                    matrix_err = true;
                }
                
                if (!matrix_err) {
                    /* Search for non-zero element to swap rows */
                    int j = q;
                    bool go = false;
                    restart = false;
                    
                    do {
                        if ((j != i) && (d[j][i] != 0.0)) {
                            /* Swap rows */
                            for (int l = 1; l <= q + 1; l++) {
                                d[i][l] = d[i][l] + d[j][l];
                            }
                            if (j > i) {
                                go = true;
                            } else {
                                restart = true;
                                i = q;
                            }
                        }
                        j--;
                        if ((!go) && (!restart) && (j == 0)) {
                            matrix_err = true;
                        }
                    } while ((j != 0) && (!go) && (!restart));
                }
                
                if (matrix_err) {
                    fatal_err("*** Matrix solution Error ***");
                }
            }
            
            if (!restart) {
                /* Normalize pivot row */
                double m = d[i][i];
                for (int l = 1; l <= q + 1; l++) {
                    d[i][l] = d[i][l] / m;
                }
                
                /* Eliminate column below and above pivot */
                for (int j = 1; j <= q; j++) {
                    if (j != i) {
                        m = d[j][i];
                        for (int l = 1; l <= q + 1; l++) {
                            d[j][l] = d[j][l] - m * d[i][l];
                        }
                    }
                }
            }
            i++;
        } while (i <= q);
    } while (restart);
    
    printf("INFO: Matrix solution completed successfully\n");
}

void output_beam(void) {
    /* Calculate and output beam forces using Pascal algorithm */
    if (!shear_flag) {
        printf("Beam    Moment-L  Moment-R  Shears-L  Shears-R\n");
        printf("(#)     (Kip-ft)  (Kip-ft)  (Kips)    (Kips)\n");
        printf("----    --------  --------  ------    ------\n");
    }
    
    for (int i = 0; i <= s - 1; i++) {
        if (!shear_flag) printf("Story %d\n", s - i);
        
        for (int j = 1; j <= p; j++) {
            int h = i * p + j;
            double o = b[h][1] / b[h][2];
            
            /* Calculate end moments using slope-deflection equations */
            double m = o * (4 * d[s + i + h][q + 1] + 2 * d[s + i + h + 1][q + 1]) - b[h][3];
            double m1 = o * (2 * d[s + i + h][q + 1] + 4 * d[s + i + h + 1][q + 1]) - b[h][4];
            
            /* Calculate shears from equilibrium */
            d[h][1] = (m + m1 + b[h][6]) / b[h][2];
            d[h][2] = (b[h][5] - m - m1) / b[h][2];
            
            if (!shear_flag) {
                printf("%11d%10.3f%10.3f%10.3f%10.3f\n", j, m, m1, d[h][1], d[h][2]);
            }
        }
    }
    if (!shear_flag) printf("\n");
}

void output_column(void) {
    /* Calculate and output column forces using Pascal algorithm */
    if (!shear_flag) {
        printf("Column   Moment    Shears    Axial\n");
        printf("(#)      (Kip-ft)  (Kips)    (Kips)\n");
        printf("----     --------  ------    ------\n");
    }
    
    for (int i = 0; i <= s - 1; i++) {
        if (!shear_flag) printf("Story %d\n", s - i);
        
        for (int j = 1; j <= p + 1; j++) {
            int h = i * (p + 1) + j;
            double o = c[h][1] / c[h][2];
            double l = 6 / c[h][2] * d[i + 1][q + 1];
            
            /* Calculate column end moments */
            double m = o * (4 * d[s + h][q + 1] + 2 * d[s + h + p + 1][q + 1] + l);
            double m1 = o * (2 * d[s + h][q + 1] + 4 * d[s + h + p + 1][q + 1] + l);
            
            if (i != s - 1) {
                l = o / c[h][2] * 6 * d[i + 2][q + 1];
                m1 = m1 - l;
                m = m - l;
            }
            
            if (j <= 1) {
                m = m - e[s - i][2];
                m1 = m1 - e[s - i][1];
            }
            
            /* Calculate axial force (accumulate from above) */
            o = 0.0;
            if (j < p + 1) o = d[i * p + j][1];
            if (j > 1) o = o + d[i * p + j - 1][2];
            if (i > 0) o = o + d[j][3];
            d[j][3] = o;
            
            /* Calculate shears */
            l = 0.0; 
            if (j == 1) l = e[s - i][3];
            
            if (!shear_flag) {
                printf("%11d%10.3f%10.3f%10.3f\n", j, m, (m + m1 - l) / c[h][2], o);
            }
            
            l = 0.0; 
            if (j == 1) l = e[s - i][4];
            
            if (!shear_flag) {
                printf("%21.3f%10.3f\n", m1, (m + m1 + l) / (-c[h][2]));
            }
        }
    }
    if (!shear_flag) printf("\n");
}

void output_rotation(void) {
    /* Output joint rotations and story displacements - exact Pascal algorithm */
    printf("Joint    Rotation     Displacement\n");
    printf("(#)      (Radians)    (Feet)\n");
    printf("----     ---------    ------\n");
    
    int l = 0;
    for (int i = s + 1; i <= q; i++) {
        int j = i - s;
        
        /* Check if at story level (displacement output) */
        if ((j - 1) % (p + 1) == 0) {  /* At story level joint */
            int story_num = s - (j - 1) / (p + 1);
            printf("Story %d\n", story_num);
            l = l + 1;
            
            if (l <= s) {
                /* Output both rotation and displacement */
                printf("%11d%15.9f%15.8f\n", j, d[i][q + 1] * 0.0048, d[l][q + 1] * 0.0048);
            } else {
                /* Output only rotation */
                printf("%11d%15.9f\n", j, d[i][q + 1] * 0.0048);
            }
        } else {
            /* Interior joint - rotation only */
            printf("%11d%15.9f\n", j, d[i][q + 1] * 0.0048);
        }
    }
    printf("\n");
}

void reset_matrix(void) {
    for (int i = 1; i <= q + 1; i++) {
        for (int j = 1; j <= q + 1; j++) {
            d[i][j] = 0.0;
        }
    }
    
    for (int i = 1; i <= MAX_MEMBERS; i++) {
        for (int j = 3; j <= 6; j++) {
            b[i][j] = 0.0;
            if (i <= MAX_WIND_STORIES) {
                e[i][j - 2] = 0.0;
            }
        }
    }
}

void user_option(void) {
    printf("============================================\n");
    printf("1. Define More Macros.\n");
    printf("2. Reload the Frame.\n");
    printf("3. Exit Program.\n");
    printf("\n");
    
    do {
        option = 0;
        printf("Enter Option ? ");
        option = readi_fixed();
        printf("%d\n", option);
    } while (option < 1 || option > 3);
    printf("\n");
}

/*=== MAIN PROGRAM ===*/
int main(void) {
    printf("<script>\n");
    printf("var last_o=null;\n");
    printf("function showTab(id_no){ let id=`tab${id_no}`; let o=document.getElementById(id); ");
    printf("o.scrollIntoView({behavior: \"auto\", block: \"nearest\", inline: \"nearest\"}); ");
    printf("o.style.backgroundColor=\"lightblue\"; ");
    printf("if (last_o != null) last_o.style.backgroundColor=\"white\"; last_o=o;}\n");
    printf("</script>\n");

    get_time_string();
    signon();

    initial_var();
    get_size();
    get_geometry();
    get_member_sizes();
    get_base_joint_type();
    option = 1;
    
    do {
        if (option == 1) {
            get_macro_loads();
        }
        apply_macro_loads();
        setup_matrix();
        solve_matrix();

        printf("\n=== ANALYSIS RESULTS ===\n");
        shear_flag = false;
        output_beam();
        output_column();
        output_rotation();

        reset_matrix();
        user_option();
    } while (option != 3);
    
    printf("*****  E N D   O F   R U N  *****\n");
    return 0;
}

/*
================================================================================
# END FILE: framer_fixed.c
#
# Key Input Parsing Fixes:
# 1. ✅ skip_comment_and_whitespace() - handles apostrophe comments properly
# 2. ✅ readr_fixed() - reads numbers while ignoring comments
# 3. ✅ readi_fixed() - reads integers while ignoring comments  
# 4. ✅ reads_fixed() - reads strings including spaces, handles comments
# 5. ✅ Proper input sequence matching Pascal exactly
# 6. ✅ Error handling that skips to next line on parse failure
#
# Expected Behavior with sample_1_story_3_bay.txt:
# - Should read "1 3" for stories/bays correctly
# - Should skip comments after apostrophes  
# - Should read beam/column sizes in correct sequence
# - Should handle load definitions properly
#
================================================================================
*/