Re: airfoil in OpenLB

[balzar29]

I’m sorry, but my c++ is very poor. I don’t understand where to change the function. I think that i’ve found the place where insert the airfoil equation, but if someone can halp me it would be better and faster. I need just some c++ explanation of this code in the file “indicatorF2D.hh”

Code:

template <typename S>
IndicatorCircle2D<S>::IndicatorCircle2D(Vector<S,2> center, S radius)
: _center(center),
_radius2(radius*radius)
{
this->_myMin = _center – radius;
this->_myMax = _center + radius;
}

// returns true if x is inside the circle
template <typename S>
bool IndicatorCircle2D<S>::operator()(bool output[], const S input[])
{
output[0] = ( std::pow(_center[0] – input[0],2) + std::pow(_center[1] – input[1], 2) <= _radius2 );
return output[0];
}

template <typename S>
bool IndicatorCircle2D<S>::distance(S& distance, const Vector<S,2>& origin, const Vector<S,2>& direction, int iC)
{
S a = direction[0]*direction[0] + direction[1]*direction[1];

// returns 0 if point is at the boundary of the sphere
if ( a == _radius2 ) {
distance = S(0);
return true;
}
// norm of direction
a = sqrt(a);

S b = 2.*((origin[0] – _center[0])*direction[0] +
(origin[1] – _center[1])*direction[1])/a;
S c = -_radius2 + (origin[0] – _center[0])*(origin[0] – _center[0])
+ (origin[1] – _center[1])*(origin[1] – _center[1]);

// discriminant
S d = b*b – 4.*c;
if (d < 0) {
return false;
}

S x1 = (- b + sqrt(d)) *0.5;
S x2 = (- b – sqrt(d)) *0.5;

// case if origin is inside the sphere
if ((x1<0.) || (x2<0.)) {
if (x1>0.) {
distance = x1;
return true;
}
if (x2>0.) {
distance = x2;
return true;
}
}
// case if origin is ouside the sphere
else {
distance = std::min(x1,x2);
return true;
}

return false;
}

template <typename S>
bool IndicatorCircle2D<S>::normal(Vector<S,2>& normal, const Vector<S,2>& origin, const Vector<S,2>& direction, int iC)
{
S dist;
if (!(distance(dist, origin, direction, iC)) ) {
return false;
}

Vector<S,2> intresection(origin + dist*direction);

normal = intresection – _center;

return true;
}

Thanks!!