1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66
|
#include <iostream>
#include <vector>
#include <cmath>
using namespace std;
//============== Geometry Stuff ========================================
struct Point{ double x, y, z; };
Point operator + ( Point p , Point q ) { return { p.x + q.x, p.y + q.y, p.z + q.z }; }
Point operator - ( Point p , Point q ) { return { p.x - q.x, p.y - q.y, p.z - q.z }; }
Point operator * ( double r, Point p ) { return { r * p.x, r * p.y, r * p.z }; }
Point operator / ( Point p , double r ) { return { p.x / r, p.y / r, p.z / r }; }
double dot ( Point p , Point q ) { return p.x * q.x + p.y * q.y + p.z * q.z ; }
Point cross ( Point p , Point q ) { return { p.y * q.z - p.z * q.y, p.z * q.x - p.x * q.z, p.x * q.y - p.y * q.x }; }
double normsq ( Point p ) { return dot( p, p ); }
double abs ( Point p ) { return sqrt( normsq( p ) ); }
ostream &operator << ( ostream &out, const Point &p ) { return out << p.x << '\t' << p.y << '\t' << p.z; }
istream &operator >> ( istream &in , Point &p ) { return in >> p.x >> p.y >> p.z; }
//======================================================================
void forceSystem( const vector<Point> &f, const vector<Point> &x, Point &F, Point &M, Point &X )
{
const double SMALL = 1.0e-30;
Point centroid;
F = M = centroid = Point{ 0, 0, 0 };
for ( int i = 0; i < f.size(); i++ )
{
F = F + f[i];
M = M + cross( x[i], f[i] );
centroid = centroid + x[i];
}
centroid = centroid / f.size();
// Find point on resultant line of action that is closest to centroid (arbitrary choice)
if ( abs( F ) < SMALL )
{
cout << "Net force is zero\n";
X = Point{ 0, 0, 0 };
}
else
{
X = ( cross( F, M ) - dot( F, centroid ) * F ) / normsq( F );
}
}
//======================================================================
int main()
{
// Data is that of your example; note the x, y, z directions and origin
vector<Point> f = { { 0, 0, -90 }, { 0, 0, -105 }, { 0, 0, -160 }, { 0, 0, -50 } };
vector<Point> x = { { 12, 1, 0 }, { 5.5, -5, 0 }, { 14.5, -3, 0 }, { 22.5, -5.5, 0 } };
Point F, M, X;
forceSystem( f, x, F, M, X );
cout << "Net force = " << F << '\n'
<< "Net moment = " << M << '\n'
<< "\nEquivalent force system is\n"
<< "force " << F << '\n'
<< "acting at " << X << '\n'
<< "(or any point on a line parallel to F passing through that point)";
}
|