Grouped weighted mean

analyze.cc

#include <cstdio>
#include <cstdlib>
#include <vector>
#include <unordered_map>
#include <utility>
using namespace std;

#include <time.h>
#include <sys/time.h>
#include <sys/resource.h>


using namespace std;

const size_t ARRAYSZ=20000000;
const size_t LINELN=100;
const size_t NBURNINS=100;
const size_t NRUNS=100;

double timespec_diff(const struct timespec&, const struct timespec&);
double timeval_diff(const struct timeval&, const struct timeval&);

int main(int argc,char *argv[])
{
    if (argc!=2){
        fprintf(stderr,"Usage: analyze [filename]\n");
        exit(1);
    }
    FILE *fp;
    fp = fopen(argv[1],"r");
    if (fp==NULL){
        perror(argv[1]);
        exit(1);
    }
    auto age=vector<char>(ARRAYSZ);
    auto sex=vector<char>(ARRAYSZ);
    auto wt=vector<int>(ARRAYSZ);

    size_t n=0;
    char line[LINELN];
    int c;
    int lineno=0;
    while(!(feof(fp)||ferror(fp))){
        if (n==ARRAYSZ){
            fprintf(stderr,"Too many lines. Skipping rest\n");
            break;
        }
        lineno++;
        int i=0;
        while ((c=fgetc(fp))!=EOF && c!='\n'){
            if (i<LINELN-1)
                line[i++] = c;
        }
        if (i==0)
            continue;
        line[i]='\0';
        int res = sscanf(line,"%hhd,%hhd,%d",&age[n],&sex[n],&wt[n]);
        if (res!=3){
            fprintf(stderr,"Error reading line no %d\n",lineno);
            fprintf(stderr,"%s\n",line);
        }
        n++;
    }
    fprintf(stderr,"%d lines read\n",lineno);

    auto utimes = vector<double>(NRUNS);
    auto stimes = vector<double>(NRUNS);
    auto ctimes = vector<double>(NRUNS);

    for (size_t r=0;r<NRUNS+NBURNINS;r++){
        struct rusage ubefore,uafter;
        struct timespec cbefore,cafter;

        if (getrusage(RUSAGE_SELF,&ubefore)){
            perror("Usage before");
            exit(1);
        }
        if (clock_gettime(CLOCK_MONOTONIC_RAW,&cbefore)){
            perror("Clock before");
            exit(1);
        }
        auto sums=unordered_map<char,pair<double,double>>();
        
        #pragma omp parallel
        {
            auto p_sums=unordered_map<char,pair<double,double>>();
            #pragma omp for
            for (size_t i=0;i<n;i++){
                auto it = p_sums.find(sex[i]);
                if (it==p_sums.end()){
                    p_sums[sex[i]] = make_pair((double)age[i]*wt[i],wt[i]);
                } else {
                    it->second.first += (double)age[i]*wt[i];
                    it->second.second += wt[i];
                }
            }
            #pragma omp critical
            {
                for (auto it=p_sums.begin();it!=p_sums.end();it++){
                    auto it2 = sums.find(it->first);
                    if (it2==sums.end()){
                        sums[it->first]=it->second;
                    } else {
                        it2->second.first += it->second.first;
                        it2->second.second += it->second.second;
                    }
                }
            }
        }
        
        if (r==0){
            for (auto it=sums.begin();it!=sums.end();it++){
                printf("%d:\t%.4f\n",it->first,
                        it->second.first/it->second.second);
            }
        }

        if (getrusage(RUSAGE_SELF,&uafter)){
            perror("Usage after");
            exit(1);
        }
        if (clock_gettime(CLOCK_MONOTONIC_RAW,&cafter)){
            perror("Clock after");
            exit(1);
        }

        if (r>=NBURNINS){
            utimes[r-NBURNINS] = timeval_diff(ubefore.ru_utime,uafter.ru_utime);
            stimes[r-NBURNINS] = timeval_diff(ubefore.ru_stime,uafter.ru_stime);
            ctimes[r-NBURNINS] = timespec_diff(cbefore,cafter);
        }
    }
    double st=0.0, ut=0.0, ct=0.0;
    for (size_t r=0;r<NRUNS;r++){
        ut += utimes[r];
        st += stimes[r];
        ct += ctimes[r];
    }
    fprintf(stderr,"Real:\t%.4fms\n",ct/NRUNS);
    fprintf(stderr,"User:\t%.4fms\n",ut/NRUNS);
    fprintf(stderr,"System:\t%.4fms\n",st/NRUNS);
}


double timeval_diff(const struct timeval& before,
                        const struct timeval& after)
{
    double btime = before.tv_sec*1e3+before.tv_usec*1.0e-3;
    double atime = after.tv_sec*1e3+after.tv_usec*1.0e-3;
    return atime-btime;
}

double timespec_diff(const struct timespec& before,
                        const struct timespec& after)
{
    double btime = before.tv_sec*1e3+before.tv_nsec*1.0e-6;
    double atime = after.tv_sec*1e3+after.tv_nsec*1.0e-6;
    return atime-btime;
}