polynomial features in C++

I use the polynomial features utility from the sklearn library in Python and now I need to convert it in C++. I found this piece of code but it does't seem to work for me :

template <class T>
std::vector<T> polynomialFeatures( const std::vector<T>& input, unsigned int degree, bool interaction_only, bool include_bias )
{
    std::vector<T> features = input;
    std::vector<T> prev_chunk = input;
    std::vector<size_t> indices( input.size() );
    std::iota( indices.begin(), indices.end(), 0 );

    for ( int d = 1 ; d < degree ; ++d )
    {
        // Create a new chunk of features for the degree d:
        std::vector<T> new_chunk;
        // Multiply each component with the products from the previous lower degree:
        for ( size_t i = 0 ; i < input.size() - ( interaction_only ? d : 0 ) ; ++i )
        {
            // Store the index where to start multiplying with the current component at the next degree up:
            size_t next_index = new_chunk.size();
            for ( auto coef_it = prev_chunk.begin() + indices[i + ( interaction_only ? 1 : 0 )] ; coef_it != prev_chunk.end() ; ++coef_it )
            {
                new_chunk.push_back( input[i]**coef_it );
            }
            indices[i] = next_index;
        }
        // Extend the feature vector with the new chunk of features:
        features.reserve( features.size() + std::distance( new_chunk.begin(), new_chunk.end() ) );
        features.insert( features.end(), new_chunk.begin(), new_chunk.end() );

        prev_chunk = new_chunk;
    }
    if ( include_bias )
        features.insert( features.begin(), 1 );

    return features;
}

I would like a simple function that is easy to use because I am a beginner in C++, I need this function in second degree, the input is 27 values and the output needs to be 405 values. Expected output for example : I call the function polyFeatures(input=[1,2,3],degree=2,interaction_only= False, include_bias= False) ==> Output=[1,2,3,1,2,3,4,6,9] ( that means input=[a,b,c] ,output=[a,b,c,a^2,ab,ac,b^2,bc,c^2]). When I execute this code in my script, I get several errors:

• first error is in std::iota( indices.begin(), indices.end(),0) : error 'iota' is not a member of 'std'.
• I also get a warning at for (auto coef_it = prev_chunk.begin()..) : warning 'auto' changes meaning in C++11, please remove it [-Wc++x-compat] and then an error : 'coef_it' does not name a type.
• an error at new_chunk.push_back(input[i]**coef_it); : invalid type argument of unary '*' (have 'int')
• and finally in main : error: no matching function to call to "polynomialFeatures(double[27],int,bool,bool)" and here because I already defined my input as double[27]