fyi benchmarks game k-nucleotide lazarus / 62-General ---------------------------------------------------------------------- igouy · 2016-05-17 23:55 The benchmarks game k-nucleotide task now explicitly requires use of a built-in / library HashMap, here's the new description - http://benchmarksgame.alioth.debian.org/u64q/knucleotide- description.html#knucleotide So this old Free Pascal program will no longer be shown - http://benchmarksgame.alioth.debian.org/u64q/program.php?test=knucleot ide&lang=fpascal&id=2 If any one is interested in contributing a program that used the library hash table implementation then here are some instructions - http://benchmarksgame.alioth.debian.org/play.html ---------------------------------------------------------------------- AlexK · 2016-05-29 18:46 I'm implementing in FreePascal same knucleotide algorithm as in C++ version by Branimir Maksimovic. I use ghashmap unit from {fpc_source}/packages/fcl-stl/src/. That library was contributed by Vlado Boza in 2011. It has one problem though. The algorithm in C++ uses custom class as a hash key, it requires overloading of "==" comparison operator for that key class to be set/get in a hash table(unique keys). This is the only comparison operator in FreePascal that can't be overloaded for classes. And ghashmap uses it, "=" for hash keys: Code: Pascal  [Select][+][-]window.onload = function(){var x1 = document.getElementById("main_content_section"); if (x1) { var x = document.getElementsByClassName("geshi");for (var i = 0; i if ((FData[h])[i]).Key=key then exit(true); For the algorithm I changed all such comparisons in ghashmap library unit to: Code: Pascal  [Select][+][-]window.onload = function(){var x1 = document.getElementById("main_content_section"); if (x1) { var x = document.getElementsByClassName("geshi");for (var i = 0; i if ((FData[h])[i]).Key.equals(key) then exit(true); The hash key custom class must provide equals method. But in this case it doesn't support simple ordinal types for keys of hash. Dilemma. ---------------------------------------------------------------------- AlexK · 2016-05-29 22:10 Even this does not work: Code: Pascal  [Select][+][-]window.onload = function(){var x1 = document.getElementById("main_content_section"); if (x1) { var x = document.getElementsByClassName("geshi");for (var i = 0; i {$mode objfpc} type   hashkey = record      data : qword; size : byte;      class operator =(lkey, rkey : hashkey) : boolean;   end; class operator hashkey.=(lkey, rkey : hashkey) : boolean;begin   result := lkey.data = rkey.data;end; begin    // Code here  end. test.pp(6,7) Fatal: Syntax error, "END" expected but "CLASS" found But this page says it should work: http://www.freepascal.org/docs- html/current/ref/refse56.html This renders FPC useless in many cases. ---------------------------------------------------------------------- molly · 2016-05-29 22:16 Sure it works: Code: [Select]program test; {$mode objfpc}{$MODESWITCH ADVANCEDRECORDS} type    hashkey = record       data : qword; size : byte;       class operator =(lkey, rkey : hashkey) : boolean;    end; class operator hashkey.=(lkey, rkey : hashkey) : boolean; begin    result := lkey.data = rkey.data; end; begin    // Code here  end. That has nothing to do with operator overloading, rather with using advanced records. As operator overloading works even without defining any mode at all (actually it defaults to mode fpc): Code: [Select]program LookNoModeDelphiOrModeFPC; type    complex = record      re : real;      im : real;    end; Var    C,Z : Complex; // New type complex    R   : Real; operator := (r : real) z : complex;  begin    z.re := r;    z.im := 0.0;  end; begin    Z := C;  // assignments between complex types.    C := R; end. You should be reading this and apply to your situation. ---------------------------------------------------------------------- AlexK · 2016-05-29 23:05 Quote from: molly on May 30, 2016, 12:16:29 am As operator overloading works even without defining any mode at all (actually it defaults to mode fpc): You should be reading this and apply to your situation. Nothing's mentioned there(official docs) about use cases for {$modeswitch advancedrecords} in a context of operator overloading. You should've pointed to a wiki page: http://wiki.freepascal.org/Record If one never used Delphi, he would've been clueless about such feature(operator overloading in advanced records). "advancedrecords" works. It helped. Quote from: molly on May 30, 2016, 12:16:29 amCode: Pascal  [Select][+][-]window.onload = function(){var x1 = document.getElementById("main_content_section"); if (x1) { var x = document.getElementsByClassName("geshi");for (var i = 0; i type    complex = record      re : real;      im : real;    end; operator := (r : real) z : complex;  begin    z.re := r;    z.im := 0.0;  end; I suppose this form of operator overloading is visible only in a unit it was defined in. Therefore it's not visible to an external library unit, and compiler says "operator not overloaded" when compiling a library generic type specialized with user's record. If a record is passed to a library, that record must contain operator overloading inside its definition. So, in this case only {$modeswitch advancedrecords} works. ---------------------------------------------------------------------- molly · 2016-05-29 23:35 QuoteNothing's mentioned there(official docs) about use cases for {$modeswitch advancedrecords} in a context of operator overloading. _You_ are the one using advanced records, ergo you need the modeswitch. Don't use advanced records and you won't need the switch. Simple as that, and as explained in the documentation about generics that you linked to. QuoteIf one never used Delphi, he would be clueless about such feature(operator overloading in advanced records). Pardon but you was the one linking to generics documentation, not advanced records. Please do not confuse the two as they are a different kind of breed. The advanced record as shown in that documentation you linked to only serves as a _type_ for which the generics class is specialized to. QuoteI suppose this form of operator overloading is visible only in a unit it was defined in. Therefore it's not visible to an external library unit, and compiler says "operator not overloaded" when compiling that library unit . Are you honestly saying that you have no idea how to declare functions and procedures in the interface section of a unit and let your other unit/program uses that unit so that unit can use those functions and procedures ? Guess what: that also works for operators ;-) QuoteIf a record is passed to a library, that record must contain operator overloading inside its definition. So only advancedrecords works in this case. Also wrong. Only if you want it to. You simply do not have to if you use classes *period* And yes, you can use class/records operators wherever you want (edit: might even be preferred) but afaik never required as you currently seem to think. edit: the reason i pointed you to generic operator overloading is that it seems the only documentation available about operator overloading, e.g. there does not seem specific documentation on using class operators and/or record operators. i was even unable to find anything on advanced records inside the online manual (but perhaps i overlooked). Thus, read what i linked to and apply for your situation. If that means you want to use record operators and/or class operators then simply use that :-) edit2: Oh, i might have been using some wrong wording there (or at least might work confusing). Yes, you need a record in case you want to enclose the operator declaration inside your class. ---------------------------------------------------------------------- AlexK · 2016-05-30 00:36 Thanks for {$modeswitch advancedrecords} hint, it's the only option that suits this algorithm. FPC does not allow overloading of equality operator for classes because it's already overloaded internally as a comparison of pointers to instances(my conclusion after reading some threads). I decided to use a record instead of a class as a hash key. Top level operator overloading is visible only in a compilation unit it was defined in. Operator overloading doesn't become a part of a record's type. If you specialize some library generic class by that record, then operator overloading for the record is ignored in methods of the specialized class. {$modeswitch advancedrecords} makes operator overloading a part of a record's type. If you specialize some library generic class by advanced record, then operator overloading for the record will be used in methods of the specialized class. Here, TSomeHashMapImpl is from some library unit, and TSomeHashMapImpl methods will use operator overloading for the hashkey record, because operator overloading of an advanced record is a part its type: Code: Pascal  [Select][+][-]window.onload = function(){var x1 = document.getElementById("main_content_section"); if (x1) { var x = document.getElementsByClassName("geshi");for (var i = 0; i {$modeswitch advancedrecords} type   hashkey = record      data : qword; size : byte;      class operator = (lkey, rkey : hashkey) : boolean;   end;   hashtable = specialize TSomeHashMapImpl; // TSomeHashMapImpl methods will use this overloadclass operator hashkey.= (lkey, rkey : hashkey) : boolean;begin   result := lkey.data = rkey.data;end; In the next code snippet TSomeHashMapImpl methods will ignore operator overloading. Error: Operator is not overloaded: "hashkey" = "hashkey": Code: Pascal  [Select][+][-]window.onload = function(){var x1 = document.getElementById("main_content_section"); if (x1) { var x = document.getElementsByClassName("geshi");for (var i = 0; i type   hashkey = record      data : qword; size : byte;   end;   hashtable = specialize TSomeHashMapImpl; // this will be used only in the current compilation unit; TSomeHashMapImpl methods do not see thisoperator = (lkey, rkey : hashkey) : boolean;begin    result := lkey.data = rkey.data;end; ---------------------------------------------------------------------- molly · 2016-05-30 01:05 QuoteHere, TSomeHashMapImpl is from some library unit, and TSomeHashMapImpl methods will use operator overloading for the hashkey record, because operator overloading of an advanced record is a part its type: Yes, that is correct AlexK. I f you want to use it the way you showed then that is perfectly ok. You are also correct about the overloading of the = or <> operators, indeed this seems impossible to do globally (the infamous "impossble operator overload" error). Ergo, in the last part that i correct you , _i_ was wrong and not you. my apologies for that.  The reason you get into problems with the last snippet you posted is because you defined a new generic using the record type. If you'd declare the hashkey and operator in the right order it should be no problem afaik. Code: [Select]type    hashkey = record       data : qword; size : byte;    end;    operator = (lkey, rkey : hashkey) : boolean; Type   hashtable = specialize TSomeHashMapImpl; Please strike that last snippet/remark. It works for me either way. main: Code: [Select]program testunits; uses   unit1; var   G1, G2: Ghashtable; begin   G1 := GhashTable.Create;   G2 := GhashTable.Create;   If G1.value = G2.value then   begin     WriteLn('are same');   end;     G1.Free;   G2.Free; end. unit1: Code: [Select]unit unit1; {$MODE OBJFPC}{$H+} interface uses   unit2; type    Thashkey = record      data : qword;      size : byte;    end; Type   Ghashtable = specialize TSomeHashMapImpl;   operator = (lkey, rkey : Thashkey) : boolean; implementation operator = (lkey, rkey : Thashkey) : boolean; begin   result := lkey.data = rkey.data; end; end. unit2: Code: [Select]unit unit2; {$MODE OBJFPC}{$H+} interface type   generic TSomeHashMapImpl = class(TObject)    public     Value : T;   end; implementation end. ---------------------------------------------------------------------- AlexK · 2016-05-30 01:39 Quote from: molly on May 30, 2016, 03:05:33 amThe reason you get into problems with the last snippet you posted is because you defined a new generic using the record type. If you'd declare the hashkey and operator in the right order it should be no problem afaik. Code: Pascal  [Select][+][-]window.onload = function(){var x1 = document.getElementById("main_content_section"); if (x1) { var x = document.getElementsByClassName("geshi");for (var i = 0; i type   hashkey = record      data : qword; size : byte;   end;    operator = (lkey, rkey : hashkey) : boolean; Type  hashtable = specialize TSomeHashMapImpl;  That won't compile. Hypothetically it could be implemented in FPC as some kind of forward declaration for an operator overload, but it's not how it should be IMO. {$modeswitch advancedrecords} is the only option. Extended record type. ---------------------------------------------------------------------- molly · 2016-05-30 01:47 it compiles perfectly fine here, and the forwarding works as advertised for me. I'm not really sure what's buggin' your code. PS: but one tends to get that with pasting half arsed snippets. The error "Error: Operator is not overloaded: "hashkey" = "hashkey":" can be solved by including your haskey definition unit in the unit with the generics declaration. So it would have been very informative telling that you are using (or at least wanted to use) the same operator inside your generics class again. ---------------------------------------------------------------------- AlexK · 2016-05-30 07:28 The difficulty of ghashmap unit consists in the fact that it does not provide default hash functor, user has to implement it. Yet it's a way to optimize hashtables for specific use cases. I just tested a primitive case statement as a hash function for small sequences(write_frequencies), and it generates the same output as a C++ version. For long sequences like "GGTATT" and "GGTATTTTAATT" a real hash function must be implemented. That's a completely different story. ---------------------------------------------------------------------- marcov · 2016-05-30 12:40 Quote from: AlexK on May 30, 2016, 01:05:51 am Nothing's mentioned there(official docs) about use cases for {$modeswitch advancedrecords} in a context of operator overloading. If you don't want to know these kind of things, simply use {$mode delphi} ---------------------------------------------------------------------- AlexK · 2016-05-30 16:11 Quote from: molly on May 30, 2016, 03:47:31 amit compiles perfectly fine here, and the forwarding works as advertised for me. You helped by pointing to {$modeswitch advancedrecords}. That worked and I moved forward. I'm sure that advancedrecords is the only way for external library units to use overloaded operators on records defined in your code, because advancedrecords contain operator overloading in their type definition, not as a top level unit definition. Quote from: molly on May 30, 2016, 03:47:31 amI'm not really sure what's buggin' your code. "Geometry".  ( https://www.youtube.com/watch?v=BKorP55Aqvg ) Quote from: molly on May 30, 2016, 03:47:31 amPS: but one tends to get that with pasting half arsed snippets. The error "Error: Operator is not overloaded: "hashkey" = "hashkey":" can be solved by including your haskey definition unit in the unit with the generics declaration. So it would have been very informative telling that you are using (or at least wanted to use) the same operator inside your generics class again. I implemented a FreePascal version of the fastest(C++) knucleotide algorithm for The Computer Language Benchmarks Game. The only problem left to solve, is to create a real hash function for ghashmap from FCL. Then multithreading maybe. I used this func in place of the real hash function temporarily: Code: Pascal  [Select][+][-]window.onload = function(){var x1 = document.getElementById("main_content_section"); if (x1) { var x = document.getElementsByClassName("geshi");for (var i = 0; i class function hashfunctor.hash(key : hashkey; b : sizeUint) : sizeUint;begin   case key.data of      0..2: result := 0;      3..5: result := 1;      6..8: result := 2;      9..11: result := 3;      12..14: result := 4;      15..255: result := 5;      256..65536: result := 6;   else      result := 7;   end;end; It returns an index of a bucket in hashmap. ---------------------------------------------------------------------- molly · 2016-05-31 02:59 First of all congratz on your progress !  Quote from: AlexK on May 30, 2016, 06:11:17 pm... I'm sure that advancedrecords is the only way for external library units to use overloaded operators on records defined in your code, because advancedrecords contain operator overloading in their type definition, not as a top level unit definition. You seem as much sure as i believe you're not entirely correct on that ;-) Since you moved one, it is of no importance anymore but, i did not wanted to let this go without showing a possibility for the casual reader. Note that i don't recommend using it, only here to show it is possible (and it might not work for your specific implementation). Of course the generic itself as well as rest of the code contains nonsense and only serves as example. Main program, includes p1unit Code: [Select]program p1testunits; uses   p1unit1; var   G1, G2  : Ghashtable;   Hashish : Array[0..1] of THashKey =   (    ( Data : 1234 ; Size: 1),    ( Data : 54321; Size: 1)   );   SamSam  : array[boolean] of string = ('differ', 'same'); begin   G1 := GhashTable.Create;   G2 := GhashTable.Create;   G1.Value := Hashish[0];   G2.Value := Hashish[0];   WriteLn('G1.Value.Data = ', G1.Value.Data);   WriteLn('G1.Value.Size = ', G1.Value.Size);   WriteLn('G2.Value.Data = ', G2.Value.Data);   WriteLn('G2.Value.Size = ', G2.Value.Size);   WriteLn('G1 and G2 ', SamSam[G1.value = G2.value]);   WriteLn;   G2.value := Hashish[1];   WriteLn('G1.Value.Data = ', G1.Value.Data);   WriteLn('G1.Value.Size = ', G1.Value.Size);   WriteLn('G2.Value.Data = ', G2.Value.Data);   WriteLn('G2.Value.Size = ', G2.Value.Size);   WriteLn('G1 and G2 ', SamSam[G1.value = G2.value]);   WriteLn;   WriteLn(G1.TestSomeThing(G1.Value));   WriteLn(G1.TestSomeThing(G2.Value));   WriteLn(G1.TestSomeThing(Hashish[0]));   WriteLn(G1.TestSomeThing(Hashish[1]));   G1.Free;   G2.Free; end. Unit 1: includes hash record, specialization of the generic and the operator. Code: [Select]unit p1unit1; {$MODE OBJFPC}{$H+} interface uses   p1unit2; type    Thashkey = record      data : qword;      size : byte;    end;   operator   = (lkey, rkey : Thashkey) : boolean; type   Ghashtable = specialize TSomeHashMapImpl; implementation operator = (lkey, rkey : Thashkey) : boolean; begin   result := lkey.data = rkey.data; end; end. Unit 2: the generic declaration. Note how unit 1 is included in the implementation section in order to use the correct operator (inside the class itself). Code: [Select]unit p1unit2; {$MODE OBJFPC}{$H+} interface type   generic TSomeHashMapImpl = class(TObject)    private     FValue : T;    public     constructor Create;     function TestSomething(SomeHash: T): string;    public     property value: T read FValue Write FValue;   end; implementation uses   p1unit1; constructor TSomeHashMapImpl.Create; begin   FValue := default(T); end; function TSomeHashMapImpl.TestSomething(SomeHash: T): string; const Fn = 'TSomeHashMapImpl.TestSomething '; begin   if (FValue = SomeHash)   then Result := Fn + 'determined that passed hash is same as internal stored hasvalue'   else Result := Fn + 'determined that passed hash differs from internal stored hasvalue'; end; end. ---------------------------------------------------------------------- AlexK · 2016-05-31 04:16 Quote from: molly on May 31, 2016, 04:59:31 amQuote from: AlexK on May 30, 2016, 06:11:17 pm... I'm sure that advancedrecords is the only way for external library units to use overloaded operators on records defined in your code, because advancedrecords contain operator overloading in their type definition, not as a top level unit definition. You seem as much sure as i believe you're not entirely correct on that ;-) ... Unit 2: the generic declaration. Note how unit 1 is included in the implementation section in order to use the correct operator (inside the class itself). That is supposed to be a library unit with a generic type. A library doesn't know about your unit. You are using that library, not library uses your unit. I explained it in initial comments. You can only specialize library generic with your type in your unit. For operator overloading to work, it must be a part of the type you specialize library generic with, not a part of a unit's interface. ---------------------------------------------------------------------- AlexK · 2016-05-31 04:29 I changed hash functor to a simple: Code: Pascal  [Select][+][-]window.onload = function(){var x1 = document.getElementById("main_content_section"); if (x1) { var x = document.getElementsByClassName("geshi");for (var i = 0; i class function hashfunctor.hash(key : hashkey; buckets : sizeUint) : sizeUint; inline;begin   result := key.data mod buckets;end; But the code is still very slow. Reading a 25 Mb file from stdin into an ansistring takes 0.6 secs on my machine. C++ does it almost immediately with std::string, 0.08 secs. C++ memory management and data structures are optimized to astounding levels. Now I'm checking if refcounting of ansistring slows things down in a function that is called many million times. Also I couldn't understand why does THashMap return 1 for the key that doesn't exist in hashmap. I tried to find something about this undefined behaviour in GetData function of THashMap: Code: Pascal  [Select][+][-]window.onload = function(){var x1 = document.getElementById("main_content_section"); if (x1) { var x = document.getElementsByClassName("geshi");for (var i = 0; i function THashmap.GetData(key:TKey):TValue;inline;var i,h,bs:longint;begin   h:=Thash.hash(key,FData.size);   bs:=(FData[h]).size;   for i:=0 to bs-1 do begin      if (((FData[h])[i]).Key=key) then exit(((FData[h])[i]).Value);   end;end; If bs is 0, then the for loop is not executed, so what is the return value of this function? What is the Result? Couldn't find in any docs for fpc or delphi. ---------------------------------------------------------------------- molly · 2016-05-31 05:02 Quote from: AlexK on May 31, 2016, 06:16:58 amThat is supposed to be a library unit with a generic type. A library doesn't know about your unit. You are using that library, not library uses your unit. Fair enough. QuoteI explained it in initial comments. It was not clear to me that the specialization was not allowed to be part of the unit containing the generic (which is another possible solution). Although i realize that would work restricting. It was not clear for me from reading the assignment descriptions alone. QuoteYou can only specialize library generic with your type in your unit. For operator overloading to work, it must be a part of the type you specialize library generic with, not a part of a unit's interface. As mentioned above, the specialization can also take part in the generics unit itself but that does work restricting in that you can't do additional specialization with using global operators (unless you incorporate those operators again). It is the nature of how globals work: you are required to include them to be able to use them. That's why i wrote that i don't recommend using it but that it is possible to do (depending on the situation). Apologies for my mix-up and the corresponding erroneous follow ups. ---------------------------------------------------------------------- AlexK · 2016-05-31 05:11 I doubt that THashMap will be able to cope with such aggressive key insertion. But in regard to the input reading, is this how it's usually done: Code: Pascal  [Select][+][-]window.onload = function(){var x1 = document.getElementById("main_content_section"); if (x1) { var x = document.getElementsByClassName("geshi");for (var i = 0; i procedure main;var   buf : shortstring;   inStr : ansistring;   //i : byte;begin   //for i := 0 to length(tochr)-1 do begin tonum[ord(tochr[i])] := i; end;   repeat readln(buf) until AnsiStartsStr('>THREE', buf); // seek   repeat      readln(buf);      if buf[0] <> ';' then inStr := inStr + buf;   until eof(input);   inStr := upcase(inStr);end; begin   main;end. From C++: Code: C  [Select][+][-]window.onload = function(){var x1 = document.getElementById("main_content_section"); if (x1) { var x = document.getElementsByClassName("geshi");for (var i = 0; i int main(){  init();  std::string input;  char buffer[256];  while (fgets(buffer,100,stdin) && memcmp(">THREE",buffer,6)!=0);  while (fgets(buffer,100,stdin) && buffer[0] != '>')    {      if (buffer[0] != ';')        {          input.append(buffer,strlen(buffer)-1);      }    }  std::transform(input.begin(),input.end(),input.begin(),::toupper);} UPD: Looked at reverse-complement benchmark, optimization is to use pchar buffers with getmem and reallocmem. ---------------------------------------------------------------------- AlexK · 2016-05-31 06:16 Quote from: molly on May 31, 2016, 07:02:58 amIt is the nature of how globals work: you are required to include them to be able to use them. That's why i wrote that i don't recommend using it but that it is possible to do (depending on the situation). Overloaded operators for records are standalone objects that can be exported in interface section, independent from the record types they operate on. Something like that. As I noted, official documentation omits this idiosyncrasy. But.. On the other hand records are not OOP, they don't receive messages, they are operated on. So this approach with operator overloading for records is explainable. {$modeswitch advancedrecords} does it the right way. Operator overloading is defined by a record's type(like classes). User can specialize imported library generic by such record and when those library generic methods are compiled in their separate unit, they "see" that the record type overrides operators. ---------------------------------------------------------------------- serbod · 2016-05-31 07:22 Try TStringHash from IniFiles unit. Memory allocation take time, allocate memory for 1000 records at once is much faster, than do 1000 allocations for each record. Do not resize arrays, strings or memory chunks - it cause data move. Use SmallString or String[nn] type, they stack-based and don't use reference counting. Use 'const' or 'var' prefix for string or record parameter, it give pointer to data instead of copy. ---------------------------------------------------------------------- serbod · 2016-05-31 12:35 My implementation from scratch: Code: Pascal  [Select][+][-]window.onload = function(){var x1 = document.getElementById("main_content_section"); if (x1) { var x = document.getElementsByClassName("geshi");for (var i = 0; i program knucleotide; uses SysUtils, Classes, IniFiles; type  TNuclCountRec = record    Name: string;    Count: Integer;  end;  TNuclCountArr = array of TNuclCountRec; var  shNucl: TStringHash;  NCArr: TNuclCountArr;  NCArrLen: Integer;  Total1, Total2: Integer; procedure CountSeq(const sSeq: String);var  n: integer;begin  n := shNucl.ValueOf(sSeq);  if n = -1 then  begin    n := NCArrLen;    Inc(NCArrLen);    SetLength(NCArr, NCArrLen);    NCArr[n].Name := sSeq;    NCArr[n].Count := 1;    shNucl.Add(sSeq, n);  end  else  begin    Inc(NCArr[n].Count);  end;end; procedure CountFrame(Len: Integer; const sFrame: ShortString);var  i, n: integer;  s: string;begin  n := Length(sFrame);  for i := 1 to Len do  begin    s := Copy(sFrame, i, 1);    CountSeq(s);    Inc(Total1);     if i < n then // last pair    begin      s := Copy(sFrame, i, 2);      CountSeq(s);      Inc(Total2);    end;     s := Copy(sFrame, i, 18);  if Pos('GGT', s) = 1 then CountSeq('GGT');    if Pos('GGTA', s) = 1 then CountSeq('GGTA');    if Pos('GGTATT', s) = 1 then CountSeq('GGTATT');    if Pos('GGTATTTTAATT', s) = 1 then CountSeq('GGTATTTTAATT');    if Pos('GGTATTTTAATTTATAGT', s) = 1 then CountSeq('GGTATTTTAATTTATAGT');  end;end; function SortItems(Item1, Item2: Pointer): Integer;begin  Result := TNuclCountRec(Item2^).Count - TNuclCountRec(Item1^).Count;end; procedure ExportFreq(n, TotalCount: Integer);var  i: Integer;  d: Double;  lst: TList;begin  lst := TList.Create();  try    for i := 0 to NCArrLen-1 do    begin      if Length(NCArr[i].Name) <> n then Continue;      lst.Add(@NCArr[i]);    end;     lst.Sort(@SortItems);     for i := 0 to lst.Count-1 do    begin      d := TNuclCountRec(lst[i]^).Count / (TotalCount / 100);    WriteLn(TNuclCountRec(lst[i]^).Name+' '+FormatFloat('0.000', d));    end;  finally    lst.Free();  end;end; procedure ExportCount(sSeq: string);var  i, n: Integer;begin  n := 0;  for i := 0 to NCArrLen-1 do  begin    if NCArr[i].Name = sSeq then    begin      n := NCArr[i].Count;      Break;    end;  end;  WriteLn(IntToStr(n) + #$09 + sSeq);end; procedure Main();var  s, sPrevFrame: ShortString;  bSkip: Boolean;  iFrameLen: Integer;begin  shNucl := TStringHash.Create(16);  NCArrLen := 0;  SetLength(NCArr, NCArrLen);  bSkip := True;  Total1 := 0;  Total2 := 0;  sPrevFrame := '';  iFrameLen := 0;  FormatSettings.DecimalSeparator := '.';   while True do  begin    ReadLn(s);    if Length(s) = 0 then Break;     if bSkip then    begin      if Pos('THREE', s) > 0 then      begin        bSkip := False;    end;      Continue;    end;     // prev frame + 18    if iFrameLen > 0 then    begin      if Length(sPrevFrame) <> (iFrameLen + 18) then      SetLength(sPrevFrame, iFrameLen + 18);      Move(s[1], sPrevFrame[iFrameLen+1], 18);      CountFrame(iFrameLen, UpperCase(sPrevFrame));    end    else    begin      SetLength(sPrevFrame, Length(s) + 18);    end;    iFrameLen := Length(s);    Move(s[1], sPrevFrame[1], iFrameLen);  end;   SetLength(sPrevFrame, iFrameLen);  CountFrame(iFrameLen, UpperCase(sPrevFrame));   shNucl.Free();   // export result  ExportFreq(1, Total1);  WriteLn();  ExportFreq(2, Total2);  WriteLn();  ExportCount('GGT');  ExportCount('GGTA');  ExportCount('GGTATT');  ExportCount('GGTATTTTAATT');  ExportCount('GGTATTTTAATTTATAGT');end; begin  Main();end. ---------------------------------------------------------------------- AlexK · 2016-05-31 18:30 Quote from: serbod on May 31, 2016, 09:22:39 amDo not resize arrays, strings or memory chunks - it cause data move. THashMap uses dynamic arrays for its key/value buckets(default memory manager). When data grows, THashMap constantly expands or contracts key/value arrays during keys rehashing. C++ libraries more likely optimize those things by direct management of memory buffers, strategy of allocation, etc. It would be interesting to see in free time how THashMap could be optimized with direct memory handling. ---------------------------------------------------------------------- AlexK · 2016-05-31 18:55 Quote from: serbod on May 31, 2016, 02:35:44 pmMy implementation from scratch: With threading(if doable) it would be close to C++ version, but still ~2 times slower. Different algorithms. On 25Mb FASTA file: C++   0m1.41s  (2 cores) FPC    0m5.38s  (1 core) ---------------------------------------------------------------------- serbod · 2016-06-01 08:55 Same with some optimizations and multithreading: Code: Pascal  [Select][+][-]window.onload = function(){var x1 = document.getElementById("main_content_section"); if (x1) { var x = document.getElementsByClassName("geshi");for (var i = 0; i { Sergey Bodrov (serbod@gmail.com) 2016-06-01 }program knucleotide; uses SysUtils, Classes, IniFiles, syncobjs; const THREADS_NUM = 4; // set max number of parallel threads type  TNuclCountRec = record    Name: string;    Count: Integer;  end;  TNuclCountArr = array of TNuclCountRec;   { TWorker }   TWorker = class(TThread)  protected    procedure Execute(); override;  public    sFrame: ShortString;    Len: Integer;    Done: Boolean;  end; var  shNucl: TStringHash;  NCArr: TNuclCountArr;  NCArrLen: Integer;  Total1, Total2: Integer;  WorkersCount: Integer;  lock: TCriticalSection;  WorkersList: TList; procedure CountSeq(const sSeq: String);var  n: integer;begin  n := shNucl.ValueOf(sSeq);  if n = -1 then  begin    lock.Acquire();    n := NCArrLen;    Inc(NCArrLen);    SetLength(NCArr, NCArrLen);    NCArr[n].Name := sSeq;    NCArr[n].Count := 1;    shNucl.Add(sSeq, n);    lock.Release();  end  else  begin    InterLockedIncrement(NCArr[n].Count);  end;end; procedure CountFrame(Len: Integer; const sFrame: ShortString);var  i, n: integer;  s: string;begin  n := Length(sFrame);  for i := 1 to Len do  begin    s := Copy(sFrame, i, 1);    CountSeq(s);    InterLockedIncrement(Total1);     if i < n then // last pair    begin      s := Copy(sFrame, i, 2);      CountSeq(s);      InterLockedIncrement(Total2);    end;     s := Copy(sFrame, i, 18);  if Pos('ggt', s) = 1 then    begin      CountSeq('ggt');      if Pos('ggta', s) = 1 then      begin        CountSeq('ggta');        if Pos('ggtatt', s) = 1 then        begin          CountSeq('ggtatt');        if Pos('ggtattttaatt', s) = 1 then          begin            CountSeq('ggtattttaatt');            if Pos('ggtattttaatttatagt', s) = 1 then CountSeq('ggtattttaatttatagt');          end;        end;      end;    end;  end;end; function SortItems(Item1, Item2: Pointer): Integer;begin  Result := TNuclCountRec(Item2^).Count - TNuclCountRec(Item1^).Count;end; procedure ExportFreq(n, TotalCount: Integer);var  i: Integer;  d: Double;  lst: TList;begin  lst := TList.Create();  try    for i := 0 to NCArrLen-1 do    begin      if Length(NCArr[i].Name) <> n then Continue;      lst.Add(@NCArr[i]);    end;     lst.Sort(@SortItems);     for i := 0 to lst.Count-1 do    begin      d := TNuclCountRec(lst[i]^).Count / (TotalCount / 100);    WriteLn(UpperCase(TNuclCountRec(lst[i]^).Name)+' '+FormatFloat('0.000', d));    end;  finally    lst.Free();  end;end; procedure ExportCount(sSeq: string);var  i, n: Integer;begin  n := 0;  for i := 0 to NCArrLen-1 do  begin    if NCArr[i].Name = sSeq then    begin      n := NCArr[i].Count;      Break;    end;  end;  WriteLn(IntToStr(n) + #$09 + UpperCase(sSeq));end; procedure SpawnWorker(Len: Integer; const sFrame: ShortString);var  Worker: TWorker;  i: Integer;begin   if WorkersList.Count < THREADS_NUM then  begin    Worker := TWorker.Create(True);    Worker.FreeOnTerminate := False;    WorkersList.Add(Worker);  end  else  begin    Worker := nil;    while WorkersCount >= THREADS_NUM do Sleep(0);    while not Assigned(Worker) do    begin      for i := 0 to WorkersList.Count-1 do      begin        if TWorker(WorkersList[i]).Done then        begin          Worker := TWorker(WorkersList[i]);        end;      end;    end;  end;   Worker.Len := Len;  Worker.sFrame := sFrame;  Worker.Done := False;  Worker.Suspended := False; end; procedure Main();var  s, sPrevFrame: ShortString;  bSkip: Boolean;  iFrameLen: Integer;  ii: Integer;begin  shNucl := TStringHash.Create(16);  lock := TCriticalSection.Create();  WorkersList := TList.Create();   NCArrLen := 0;  SetLength(NCArr, NCArrLen);  bSkip := True;  Total1 := 0;  Total2 := 0;  sPrevFrame := '';  iFrameLen := 0;  FormatSettings.DecimalSeparator := '.';  WorkersCount := 0;  ii := 0;   while True do  begin    ReadLn(s);    if Length(s) = 0 then Break;     if bSkip then    begin      if Pos('>THREE', s) > 0 then    begin        bSkip := False;      end;      Continue;    end;     Inc(ii);     // prev frame + 18    if iFrameLen > 0 then    begin    if Length(sPrevFrame) <> (iFrameLen + 18) then        SetLength(sPrevFrame, iFrameLen + 18);      Move(s[1], sPrevFrame[iFrameLen+1], 18);      if THREADS_NUM > 1 then        SpawnWorker(iFrameLen, sPrevFrame)      else        CountFrame(iFrameLen, sPrevFrame);    end    else    begin      SetLength(sPrevFrame, Length(s) + 18);    end;    iFrameLen := Length(s);    Move(s[1], sPrevFrame[1], iFrameLen);  end;   SetLength(sPrevFrame, iFrameLen);  CountFrame(iFrameLen, sPrevFrame);   while WorkersCount > 0 do Sleep(1);   WorkersList.Free();  lock.Free();  shNucl.Free();   // export result  ExportFreq(1, Total1);  WriteLn();  ExportFreq(2, Total2);  WriteLn();  ExportCount('ggt');  ExportCount('ggta');  ExportCount('ggtatt');  ExportCount('ggtattttaatt');  ExportCount('ggtattttaatttatagt');end; { TWorker } procedure TWorker.Execute;begin  while not Terminated do  begin    if not Done then    begin      InterlockedIncrement(WorkersCount);      CountFrame(Len, sFrame);      InterlockedDecrement(WorkersCount);    Done := True;    end;    Sleep(0);  end;end; begin  Main();end. ---------------------------------------------------------------------- AlexK · 2016-06-01 20:43 Quote from: serbod on June 01, 2016, 10:55:48 amSame with some optimizations and multithreading It's 2.1 times slower than your first sequential version, for some reason. I added cthreads unit in the uses clause to compile it. But faster than my port of the C++ version. It gives slightly distorted results for long strings due to a critical bug in ghashmap unit, that I fixed locally(one line). It's slow in memory management, tens of thousands of dynamic arrays are constantly resized in library implementation. Would be interesting to optimize ghashmap some time. Code: Pascal  [Select][+][-]window.onload = function(){var x1 = document.getElementById("main_content_section"); if (x1) { var x = document.getElementsByClassName("geshi");for (var i = 0; i {Alex Karbivnichiy} {$mode objfpc}{$modeswitch advancedrecords} uses sysutils, strutils, ghashmap, gmap, gutil; const tochr : array[0..3] of char = ('A', 'C', 'T', 'G');var   tonum : array[0..255] of byte; type   hashkey = record      data : qword; size : byte;      class operator = (lkey, rkey : hashkey) : boolean;   end;   hashfunctor = class   public class function hash(key : hashkey; buckets : sizeUint) : sizeUint; inline;   end;   hashtable = specialize THashMap; class operator hashkey.= (lkey, rkey : hashkey) : boolean;begin   result := lkey.data = rkey.data;end; procedure reset(var hk : hashkey; inStr : ansistring; beg : longword; fin : longword); inline;begin   hk.size := fin - beg;   hk.data := 0;   repeat      hk.data := hk.data shl 2;      hk.data := hk.data or qword(tonum[ord(inStr[beg+1])]);      inc(beg);   until beg = fin;end; function keyToStr(hk : hashkey) : ansistring;var   tmp : ansistring = '';   tmp1 : qword;   i : longword = 0;begin   tmp1 := hk.data;   while i <> hk.size do begin      tmp := tmp + tochr[tmp1 and 3];      tmp1 := tmp1 shr 2;      inc(i);   end;   result := tmp;end; class function hashfunctor.hash(key : hashkey; buckets : sizeUint) : sizeUint; inline;begin   result := key.data mod buckets;end; function calculate(inStr : ansistring; size : byte; beg : byte = 0) : hashtable;var   frequencies : hashtable;   k : hashkey;   i : integer;begin   frequencies := hashtable.create;   for i:=beg to length(inStr)+1-size do begin      reset(k, inStr, i, i+size);      frequencies[k] := frequencies[k] + 1;   end;   result := frequencies;end; procedure writeFrequencies(inStr : ansistring; size : byte);type   stgreater = specialize TGreater;   stmap = specialize TMap;var   sum : longword;   frequencies: hashtable;   freq : stmap;   iter1 : hashtable.TIterator;   iter2 : stmap.TIterator;   psum : double;begin   sum := length(inStr) + 1 - size;   frequencies := calculate(inStr, size);   freq := stmap.create;   iter1 := frequencies.iterator;   repeat      freq.insert(iter1.value, keyToStr(iter1.key));   until not iter1.next;   iter2 := freq.min;   repeat      if sum <> 0 then psum := double(100 * iter2.key) / sum else psum := 0.0;      writeln(format('%s  %-2.3f', [iter2.value, psum]));   until not iter2.next;   writeln;end; procedure writeCount(inStr : ansistring; str : ansistring);var   size : byte;   frequencies : hashtable;   k : hashkey;begin   size := length(str);   frequencies := calculate(inStr, size);   reset(k, str, 0, size);    writeln(format('%-7d  %s', [frequencies[k], str]));end; procedure main;var   buf : shortstring;   inStr : ansistring;   i : byte;begin   for i := 0 to length(tochr)-1 do begin tonum[ord(tochr[i])] := i; end;   repeat readln(buf) until AnsiStartsStr('>THREE', buf); // seek   repeat      readln(buf);      if buf[1] <> ';' then inStr := inStr + buf;   until eof(input);      inStr := upcase(inStr);    writeFrequencies(inStr, 1);   writeFrequencies(inStr, 2);   writeCount(inStr, 'GGT');   writeCount(inStr, 'GGTA');   writeCount(inStr, 'GGTATT');   writeCount(inStr, 'GGTATTTTAATT');   writeCount(inStr, 'GGTATTTTAATTTATAGT');end; begin   main;end.  Critical bug in ghashmap, added line 9: Code: Pascal  [Select][+][-]window.onload = function(){var x1 = document.getElementById("main_content_section"); if (x1) { var x = document.getElementsByClassName("geshi");for (var i = 0; i function THashmap.GetData(key:TKey):TValue;inline;var i,h,bs:longint;begin   h:=Thash.hash(key,FData.size);   bs:=(FData[h]).size;   for i:=0 to bs-1 do begin      if (((FData[h])[i]).Key=key) then exit(((FData[h])[i]).Value);   end;   result := default(TValue); // added this lineend; ---------------------------------------------------------------------- AlexK · 2016-06-07 03:32 I rewrote THashMap, and now the test shows same performance as C++. I still have to refactor some code. One interesting discovery, -gv option gives additional ~11% speedup to the -O4 option. UPD: -gv enforces cmem unit, so that's why. ---------------------------------------------------------------------- serbod · 2016-06-07 08:39 QuoteThe work is to use the ordinary built-in or library hash table implementation to accumulate count values - retrieve the count value for a key, update the count value for that key. Don't optimize away the work. ---------------------------------------------------------------------- Leledumbo · 2016-06-07 09:14 Quote from: serbod on June 07, 2016, 10:39:43 amQuoteThe work is to use the ordinary built-in or library hash table implementation to accumulate count values - retrieve the count value for a key, update the count value for that key. Don't optimize away the work. Doesn't say that optimizing the underlying library is not allowed If AlexK submit his THashMap optimization as patches, then we will get the optimization for free. ---------------------------------------------------------------------- AlexK · 2016-06-20 17:11 C++ version virtually creates a binary Trie(https://en.wikipedia.org/wiki/Trie). Uses an optimized collisionless hash function for this particular use case(only 4 letters), and uses GCC's hashmap library as an expandable array(FPC has dynamic arrays for that). I switched to creating a more general character Trie in FPC, will come back to it later.